Mid-Cayman Rise biomass analysis
Sarah Hu
2/3/2021
RMarkdown history:
Feb 2021 started github repo (https://github.com/shu251/midcayman-rise-microeuk)
continue to commit changes and as you include knitr and commit changes, those will be pushed to the git repo
eventually set up the gitpages under settings
https://resources.github.com/whitepapers/github-and-rstudio/
1 Mid-Cayman Rise grazing
Summary of cruise operations
Below is the compilation of all microscopy data
1.1 Set up working environment
GitHub repo hosts all raw data to reproduce analysis and results. See input-data or clone the repo.
library(tidyverse); library(cowplot); library(broom)2 Protistan cell concentration
Import eukaryotic cell count data from grazing experiments. In this section, we will calculate cells per ml from raw counts (Field of view, etc.) and use to estimate protist cell concentration. These will be used below in grazing experiment calculations.
counts <- read.delim("input-data/euk-counts-compiled.txt",
blank.lines.skip = FALSE,
na.strings = c("", "NA"),
stringsAsFactors = FALSE) # Import
counts[is.na(counts)] <- 0 # Change blanks to zeroesRaw data table collected during microscopy count process. Below code reviews the structure of this raw data and updates column headers to be more ‘R’ friendly.
colnames(counts) <- c("DATE", "SAMPLE", "EXPID", "VOL", "MAG", "FOV", "nanoNoFLP", "microNoFLP", "nanoFLP", "microFLP", "NOTES", "DateCompiled"); colnames(counts)## [1] "DATE" "SAMPLE" "EXPID" "VOL" "MAG"
## [6] "FOV" "nanoNoFLP" "microNoFLP" "nanoFLP" "microFLP"
## [11] "NOTES" "DateCompiled"# Column 9 and 10 list the total number of FLP inside each, so # of occurences
# Notes if count was not countable
# unique(counts$NOTES)
head(counts)## DATE SAMPLE EXPID VOL MAG FOV nanoNoFLP microNoFLP nanoFLP microFLP
## 1 4/2/20 VD-Plume T0-Rep1 150 100 1 0 0 0 0
## 2 4/2/20 VD-Plume T0-Rep1 150 100 2 0 0 0 0
## 3 4/2/20 VD-Plume T0-Rep1 150 100 3 3 0 0 0
## 4 4/2/20 VD-Plume T0-Rep1 150 100 4 1 0 4 0
## 5 4/2/20 VD-Plume T0-Rep1 150 100 5 0 1 0 0
## 6 4/2/20 VD-Plume T0-Rep1 150 100 6 0 0 0 0
## NOTES DateCompiled
## 1 0 6/7/20
## 2 0 6/7/20
## 3 0 6/7/20
## 4 0 6/7/20
## 5 elongated cell? 6/7/20
## 6 0 6/7/20To count occurence and number of FLP ingested by eukaryotic cells, the number of FLPs ingested was tallied and comma separated for multiple eukaryotic cells with FLP. These values need to separated and counted as 1 eukaryotic cell each, but retain the number of FLP per cell.
counts_occur <- counts %>%
# remove incomplete
filter(NOTES != "Not countable") %>%
# Count number of euk cells observed with FLPs (ex. if "1,2", 'occur' will = 2)
mutate(nanoFLP_occur = as.numeric(str_count(nanoFLP, "[1-9]\\d*")),
microFLP_occur = as.numeric(str_count(microFLP, "[1-9]\\d*")),
# Add number of euk cells with FLPs to those without for total number of euk cells
nanoTOTAL = as.numeric(nanoNoFLP) + nanoFLP_occur,
microTOTAL = as.numeric(microNoFLP) + microFLP_occur,
euksTOTAL = nanoTOTAL + microTOTAL) %>%
data.frame
# unique(counts$nanoFLP)
# unique(counts_occur$nanoFLP_occur)2.1 Calculate cells per ml (euk)
Input data are the raw microscopy counts by FOV. Code below calculations cells/ml based on these values. Additionaly, variance and standard deviation are also calculated. Eukaryotic cells were also classified by size, where micro equates to >20um and nano is <20um. All counts were done at 100x magnification, confirm this:
unique(counts_occur$MAG)## [1] 100# head(counts_occur)Calculate cell concentration (cells/ml).
counts_cellsml_all <- counts_occur %>%
group_by(SAMPLE, EXPID, VOL) %>% #Calculate averages by sample
summarise(totalFOV = n(), # Count total FOV counted
nanoAvg = sum(nanoTOTAL)/totalFOV, #Average per FOV
nanoVar = var(nanoTOTAL), #Variance
nanoSd = (2*(sqrt(nanoVar))), #Standard deviation
microAvg = sum(microTOTAL)/totalFOV, ## Repeat for microeuks
microVar = var(microTOTAL),
microSd = (2*(sqrt(microVar))),
euksAvg = sum(euksTOTAL)/totalFOV, ## Repeat for total cell count
euksVar = var(euksTOTAL),
euksSd = (2*(sqrt(euksVar))),
.groups = 'drop_last') %>%
# Calculate cells/ml based on magnification (at x100, 0.01 is vol of grid), volume filtered (VOL), dilution factor (0.9), and area of counting grid (for Huber lab scope, it is 283.385):
mutate(nanoCONC = ((nanoAvg * 283.385)/(VOL * 0.01 * 0.9)),
microCONC = ((microAvg * 283.385)/(VOL * 0.01 * 0.9)),
eukCONC = ((euksAvg * 283.385)/(VOL * 0.01 * 0.9))
) %>%
# left_join(expmeta) %>%
separate(SAMPLE, c("Site", "Name"), sep = "-", remove = FALSE) %>%
separate(EXPID, c("TimePoint", "Replicate"), sep = "-", remove = FALSE) %>%
data.framehead(counts_cellsml_all)## SAMPLE Site Name EXPID TimePoint Replicate VOL
## 1 Piccard-LotsOShrimp Piccard LotsOShrimp T0-Rep3 T0 Rep3 50
## 2 Piccard-LotsOShrimp Piccard LotsOShrimp T15-Rep3 T15 Rep3 50
## 3 Piccard-LotsOShrimp Piccard LotsOShrimp T20-Rep3 T20 Rep3 50
## 4 Piccard-LotsOShrimp Piccard LotsOShrimp T40-Rep3 T40 Rep3 40
## 5 Piccard-Plume Piccard Plume T0-Rep1 T0 Rep1 150
## 6 Piccard-Plume Piccard Plume T0-Rep2 T0 Rep2 150
## totalFOV nanoAvg nanoVar nanoSd microAvg microVar microSd
## 1 30 0.3666667 0.3781609 1.2298958 0.00000000 0.00000000 0.0000000
## 2 30 0.3000000 0.2172414 0.9321832 0.06666667 0.06436782 0.5074163
## 3 30 0.3333333 0.2988506 1.0933445 0.00000000 0.00000000 0.0000000
## 4 30 0.1666667 0.1436782 0.7580980 0.00000000 0.00000000 0.0000000
## 5 30 0.3000000 0.2862069 1.0699662 0.03333333 0.03333333 0.3651484
## 6 30 0.2000000 0.1655172 0.8136762 0.06666667 0.06436782 0.5074163
## euksAvg euksVar euksSd nanoCONC microCONC eukCONC
## 1 0.3666667 0.3781609 1.229896 230.90630 0.00000 230.90630
## 2 0.3666667 0.2402299 0.980265 188.92333 41.98296 230.90630
## 3 0.3333333 0.2988506 1.093345 209.91481 0.00000 209.91481
## 4 0.1666667 0.1436782 0.758098 131.19676 0.00000 131.19676
## 5 0.3333333 0.3678161 1.212957 62.97444 6.99716 69.97160
## 6 0.2666667 0.2712644 1.041661 41.98296 13.99432 55.97728Replicates belong to the same experiment for either Bag or IGT incubation. Below, modify these names and label new column with bag or igt. And create an average across replicates.
Average cells/ml across replicates, pivot to long format
counts_cellsml_avg <- counts_cellsml_all %>%
select(Site, Name, TimePoint, Replicate, nanoCONC, microCONC, eukCONC) %>%
mutate(EXP_TYPE = case_when(
grepl("IGT", Replicate) ~ "IGT",
grepl("Rep", Replicate) ~ "Bag"
)) %>%
mutate(IGT_REP = case_when(
EXP_TYPE == "IGT" ~ Replicate,
EXP_TYPE == "Bag" ~ "Bag")) %>%
select(-Replicate) %>%
pivot_longer(cols = ends_with("CONC"), names_to = "VARIABLE", values_to = "CONCENTRATION") %>%
group_by(Site, Name, TimePoint, EXP_TYPE, IGT_REP, VARIABLE) %>%
# Calculate mean, variance, SD, min, and max
summarise(MEAN = mean(CONCENTRATION),
VAR = var(CONCENTRATION),
SD = sd(CONCENTRATION),
SEM =(sd(CONCENTRATION)/sqrt(length(CONCENTRATION))),
MIN = min(CONCENTRATION),
MAX = max(CONCENTRATION),
.groups = 'drop_last') %>%
data.frame
head(counts_cellsml_avg)## Site Name TimePoint EXP_TYPE IGT_REP VARIABLE MEAN VAR SD SEM
## 1 Piccard LotsOShrimp T0 Bag Bag eukCONC 230.90630 NA NA NA
## 2 Piccard LotsOShrimp T0 Bag Bag microCONC 0.00000 NA NA NA
## 3 Piccard LotsOShrimp T0 Bag Bag nanoCONC 230.90630 NA NA NA
## 4 Piccard LotsOShrimp T15 Bag Bag eukCONC 230.90630 NA NA NA
## 5 Piccard LotsOShrimp T15 Bag Bag microCONC 41.98296 NA NA NA
## 6 Piccard LotsOShrimp T15 Bag Bag nanoCONC 188.92333 NA NA NA
## MIN MAX
## 1 230.90630 230.90630
## 2 0.00000 0.00000
## 3 230.90630 230.90630
## 4 230.90630 230.90630
## 5 41.98296 41.98296
## 6 188.92333 188.92333# View(counts_cellsml_avg)
# NOTES on calculations:
# VAR = takes the sum of the squares of each value's deviation from the mean and divides by the number of such values minus one. This differs from the calculation of variance across an entire population in that the latter divides by the size of the dataset without subtracting one.
# SD = standard deviation of all values
# SEM = standard deviation of sampling distribution; standard deviation divided by the square root of the sample size.# Euk cell counts, all and averaged
# head(counts_cellsml_all)
# head(counts_cellsml_avg)Parse experiment type information
# Convert to long format and add column that reports IGT vs bag experiment
plot_euk_conc <- counts_cellsml_all %>%
select(Site, Name, TimePoint, Replicate, ends_with("CONC")) %>%
mutate(EXP_TYPE = case_when(
grepl("IGT", Replicate) ~ "IGT",
grepl("Rep", Replicate) ~ "Bag"
)) %>%
pivot_longer(cols = ends_with("CONC"), names_to = "VARIABLE", values_to = "CONCENTRATION") %>%
data.frame# head(plot_euk_conc)
unique(plot_euk_conc$Name)## [1] "LotsOShrimp" "Plume" "Shrimpocalypse" "BSW"
## [5] "MustardStand" "OMT" "Rav2" "ShrimpHole"
## [9] "X18"vent_ids <- c("BSW","Plume", "Shrimpocalypse", "LotsOShrimp", "X18", "OMT", "Rav2", "MustardStand", "ShrimpHole")
vent_fullname <- c("Background","Plume", "Shrimpocalypse", "Lots 'O Shrimp", "X-18", "Old Man Tree", "Ravelin #2", "Mustard Stand", "Shrimp Hole")
site_ids <- c("VD", "Piccard")
site_fullname <- c("Von Damm", "Piccard")
plot_euk_conc$SiteOrder <- factor(plot_euk_conc$Site, levels = site_ids, labels = site_fullname)
plot_euk_conc$NameOrder <- factor(plot_euk_conc$Name, levels = vent_ids, labels = vent_fullname)2.2 Plot eukaryote cells per ml
Box plot to report all eukaryote cell counts.
# Code for plot
conc_boxplot <- ggplot(plot_euk_conc, aes(x = NameOrder,
y = CONCENTRATION,
group = NameOrder,
fill = VARIABLE,
shape = EXP_TYPE)) +
geom_boxplot() +
# Do not color by time point
geom_jitter(color = "black", size = 2, aes(fill = VARIABLE,
shape = EXP_TYPE)) +
scale_shape_manual(values = c(21,24)) +
scale_fill_manual(values = c("#e7298a", "#fcbba1", "#c6dbef")) +
coord_flip() +
scale_y_log10() +
# scale_y_log10(limits = c(10,1000), expand = c(0, 0)) +
facet_grid(SiteOrder ~ EXP_TYPE, space = "free", scale = "free") +
theme_bw() +
theme(axis.text.x = element_text(angle = 0, h = 1, vjust = 1),
strip.background = element_blank(),
legend.position = "right",
legend.title = element_blank()) +
labs(x = "", y = bquote("Eukaryote cells "~mL^-1),
title = "Distribution of all eukaryotic cell counts")Eukaryote cell concentration (cells/ml) are lower in the background and plume samples compared to vent sites. ~300 cells/ml in background and plume compared to ~1000 cells per ml at the vent sites. These values are also consistent between each vent site (Von Damm and Piccard) and between Bag and IGT samples.
Boxplot represents the median (line in box) and the 1st and 3rd quartiles in the lower and upper hinges, respectively (25th and 75th percentiles). Black data points are outliers from the boxplot. Upper and lower whiskers represent the 1.5 * interquartile ranges. Pink data points are the values contributing to the boxplot (individial counts across replicates and time points.)
conc_boxplot## Warning: Transformation introduced infinite values in continuous y-axis
## Warning: Transformation introduced infinite values in continuous y-axis## Warning: Removed 39 rows containing non-finite values (stat_boxplot).## Warning: Removed 39 rows containing missing values (geom_point).
eukCONC is the sum of micro and nano. Because there was a discrepency between the micro and nano cell counts, we plan to combine for most of the analysis. Here we show that the cell concentration across replicate samples was similar throughout experiments. And that the bag versus IGT experiment results were within range of one another.
Subset to plot from T0 only.
vent_ids <- c("BSW","Plume", "Shrimpocalypse", "LotsOShrimp", "X18", "OMT", "Rav2", "MustardStand", "ShrimpHole")
vent_fullname <- c("Background","Plume", "Shrimpocalypse", "Lots 'O Shrimp", "X-18", "Old Man Tree", "Ravelin #2", "Mustard Stand", "Shrimp Hole")
plot_euk_format <- plot_euk_conc %>%
filter(TimePoint == "T0" & (VARIABLE == "eukCONC")) %>%
group_by(SiteOrder, NameOrder, TimePoint, EXP_TYPE, VARIABLE) %>%
summarise(avg_conc = mean(CONCENTRATION),
SEM_conc = (sd(CONCENTRATION)/sqrt(length(CONCENTRATION))),
.groups = "rowwise") %>%
unite(EXPERIMENT, SiteOrder, NameOrder, EXP_TYPE, remove = FALSE) %>%
data.frame
# Factor
plot_euk_format$Site_Order <- factor(plot_euk_format$SiteOrder, levels = site_fullname, labels = site_fullname)
# View(plot_euk_format)
euk_plot <- ggplot(plot_euk_format, aes(x = NameOrder, y = avg_conc, fill = Site_Order)) +
geom_errorbar(aes(ymax = (avg_conc + SEM_conc), ymin = (avg_conc - SEM_conc)), width = 0.2) +
geom_point(aes(fill = Site_Order), color = "black", stat = "identity", size = 3, shape = 23) +
facet_grid(.~ Site_Order, space = "free", scales = "free") +
scale_fill_manual(values = c("#1c9099", "#de2d26")) +
theme_minimal() +
theme(panel.grid.major = element_line(), panel.grid.minor = element_blank(),
panel.background = element_blank(),
axis.line = element_line(colour = "black"),
axis.text.x = element_text(color="black", size = 12,
angle = 45, hjust = 1, vjust = 1),
axis.text.y = element_text(color="black", size = 12),
axis.title =element_text(color="black", size = 12),
axis.ticks = element_line(),
strip.text =element_blank(), legend.title = element_blank()) +
labs(x = "", y = bquote("Eukaryote cells "~mL^-1),
title = "")euk_plot
Summary of findings from eukaryote cell concentration.
# head(plot_euk_format)
plot_euk_format %>%
type.convert(as.is = TRUE) %>%
filter(VARIABLE == "eukCONC") %>%
mutate(SAMPLE_TYPE = case_when(
NameOrder == "Background" ~ "Background",
NameOrder == "Plume" ~ "Plume",
TRUE ~ SiteOrder
)) %>%
group_by(SAMPLE_TYPE) %>%
summarise(MEAN_cellml = format(mean(avg_conc), scientific = T),
min_cellml = format(min(avg_conc), scientific = T),
max_cellml = format(max(avg_conc), scientific = T),
num = n())## # A tibble: 4 × 5
## SAMPLE_TYPE MEAN_cellml min_cellml max_cellml num
## <chr> <chr> <chr> <chr> <int>
## 1 Background 9.183773e+01 9.183773e+01 9.183773e+01 1
## 2 Piccard 3.801791e+02 2.309063e+02 4.548154e+02 3
## 3 Plume 1.185379e+02 7.930115e+01 1.577747e+02 2
## 4 Von Damm 4.105001e+02 2.597696e+02 6.20998e+02 63 Prokaryote cell concentration
Import counts from DAPI counts of bacteria and archaea.
prok <- read.delim("input-data/prokINSITU-counts-compiled.txt")
# head(prok)
# unique(prok$CELLML)
insitu_proks <- prok %>%
filter(CELLML != "not countable") %>%
separate(SAMPLE, c("Site", "Name"), sep = "-", remove = FALSE) %>%
group_by(SAMPLE, Site, Name) %>%
summarise(MEAN = mean(as.numeric(CELLML)),
SD = sd(CELLML),
SEM = (sd(CELLML)/sqrt(length(CELLML))),
.groups = "rowwise") %>%
data.frame
# head(insitu_proks)
# View(insitu_proks)3.1 Microscopy cells per ml (prok)
Plot in situ prokaryote values. Average across dives where applicable.
insitu_proks$Name_order <- factor(insitu_proks$Name, levels = c("BSW", "Plume", "Quakeplume", "Shrimpocalypse", "LotsOShrimp", "X18", "OMT", "Rav2", "MustardStand", "ShrimpHole", "HotChimlet1", "ShrimpGulley", "SouthofHotChimlet", "SouthofLungSnack", "ArrowLoop", "Bartizan", "Rav1"), labels = c("Background","Plume", "Quakeplume", "Shrimpocalypse", "Lots 'O Shrimp", "X-18", "Old Man Tree", "Ravelin #2", "Mustard Stand", "Shrimp Hole", "Hot Chimlet #1", "Shrimp Gulley", "South of Hot Chimlet", "South of LungSnack", "Arrow Loop", "Bartizan", "Ravelin #1"))
site_ids <- c("VD", "Piccard")
site_fullname <- c("Von Damm", "Piccard")
insitu_proks$Site_order <- factor(insitu_proks$Site, levels = site_ids, labels = site_fullname)
# unique(insitu_proks$Name)
prok_plot <- ggplot(insitu_proks, aes(x = Name_order, y = MEAN)) +
geom_errorbar(aes(ymax = (MEAN + SEM), ymin = (MEAN - SEM)), width = 0.2) +
geom_point(stat = "identity", shape = 23, aes(fill = Site), size = 3) +
facet_grid(.~ Site_order, space = "free", scales = "free") +
scale_fill_manual(values = c("#de2d26", "#1c9099")) +
labs(y = bquote("Prokaryote cells "~mL^-1), x = "", title = "") +
scale_y_log10() +
theme_minimal() +
theme(panel.grid.major = element_line(), panel.grid.minor = element_blank(),
panel.background = element_blank(),
axis.line = element_line(colour = "black"),
axis.text.x = element_text(color="black", size = 12,
angle = 45, hjust = 1, vjust = 1),
axis.text.y = element_text(color="black", size = 12),
axis.title =element_text(color="black", size = 12),
axis.ticks = element_line(),
strip.text =element_blank(), legend.title = element_blank())prok_plot
Summary of findings from prokaryote cell counts.
insitu_proks %>%
filter(Name != "Quakeplume") %>%
mutate(SAMPLE_TYPE = case_when(
Name == "BSW" ~ "Background",
Name == "Plume" ~ "Plume",
TRUE ~ Site
)) %>%
group_by(SAMPLE_TYPE) %>%
summarise(MEAN_cellml = format(mean(MEAN), scientific = T),
min_cellml = format(min(MEAN), scientific = T),
max_cellml = format(max(MEAN), scientific = T),
num = n())## # A tibble: 4 × 5
## SAMPLE_TYPE MEAN_cellml min_cellml max_cellml num
## <chr> <chr> <chr> <chr> <int>
## 1 Background 2.487491e+04 1.186019e+04 3.788962e+04 2
## 2 Piccard 1.097126e+05 5.387814e+04 2.385857e+05 6
## 3 Plume 3.395372e+04 1.647831e+04 5.142913e+04 2
## 4 VD 4.090657e+04 8.816422e+03 1.114298e+05 63.2 Comparison with previous MCR
Compare in situ prokaryote cell counts from 2020 to previous years
prok_prev <- read.csv("input-data/cellcount_previousyr.csv")
# View(prok_prev)
# head(prok_prev)
prok_prev_formatted <- prok_prev %>%
mutate(VENTSITE = case_when(
grepl("Piccard", Site) ~ "Piccard",
grepl("Von Damm", Site) ~ "VD"
)) %>%
filter(!is.na(YEAR)) %>% #QC of
# filter(cells_ml != "NC") %>%
# filter(cells_ml != "") %>%
# filter(cells_ml != "no data") %>%
type.convert(as.is = TRUE, numerals = "no.loss") %>%
select(YEAR, VENTSITE, NAME = Name, REP=Replicate, CELLML = cells_ml, ORIGSAMPLE = Orig_vent_site_ID, ID_number, Origin)
# str(prok_prev_formatted)
# View(prok_prev_formatted)
# unique(prok_prev_formatted$VENTSITE)# Re-import 2020
prok <- read.delim("input-data/prokINSITU-counts-compiled.txt")
# View(prok)
proks_allyrs <- prok %>%
separate(SAMPLE, c("VENTSITE", "NAME"), sep = "-", remove = FALSE) %>%
mutate(YEAR = 2020) %>%
select(YEAR, VENTSITE, NAME, REP, CELLML, ORIGSAMPLE = BAC) %>%
bind_rows(prok_prev_formatted %>% select(-ID_number, -Origin)) %>%
type.convert(as.is = TRUE) %>%
# Remove not countable or not data samples:
filter(CELLML != "NC") %>%
filter(CELLML != "") %>%
filter(CELLML != "no data") %>%
filter(CELLML != "not countable") %>%
data.frame
# View(proks_allyrs)
# View(as.data.frame(unique(proks_allyrs$NAME)))
vent_order <- c("BSW","Plume","Quakeplume","NearsummitBeebee","MainOrifice","NearMainOrifice","Rav1","HotChimlet1","HotChimlet","SouthofHotChimlet","NearHotChimlet","HotCracks1","HotCracks2","ShrimpHole","ShrimpHole(X18)","X18","X19","SouthofLungSnack","TwinPeaks","OMT","WhiteCastle","GingerCastle","ArrowLoop","Bartizan","LotsOShrimp","MustardStand","ShrimpButtery","ShrimpCanyon","ShrimpGulley","Shrimpocalypse","ShrimpVegas")
vent_names <- c("Background","Plume","Quakeplume","Near summit Beebee Vents Mound","Main Orifice","Near Main Orifice","Ravelin #1","Hot Chimlet #1","Hot Chimlet","South of Hot Chimlet","Near Hot Chimlet","Hot Cracks #1","Hot Cracks #2","Shrimp Hole","Shrimp Hole (X-18)","X-18","X-19","South of Lung Snack","Twin Peaks","Old Man Tree","White Castle","Ginger Castle","Arrow Loop","Bartizan","Lots O Shrimp","Mustard Stand","Shrimp Buttery","Shrimp Canyon","Shrimp Gulley","Shrimpocalypse","Shrimp Vegas")
proks_allyrs$NAME_ORDER <- factor(proks_allyrs$NAME, levels = vent_order, labels = vent_names)
proks_allyrs$VENTSITE_ORDER <- factor(proks_allyrs$VENTSITE, levels = c("Piccard", "VD"), labels = c("Piccard", "Von Damm"))# pdf("compare-across-yr-cellcount-04052021.pdf", h = 8, w = 7)
ggplot(proks_allyrs, aes(x = NAME_ORDER, y = as.numeric(CELLML), fill = factor(YEAR), shape = VENTSITE_ORDER)) +
geom_point(stat = "identity", aes(fill = factor(YEAR)), size = 3) +
scale_shape_manual(values = c(21,23)) +
coord_flip() +
facet_grid(VENTSITE_ORDER ~ ., space = "free", scales = "free") +
scale_y_log10() +
scale_fill_manual(values = c("#1c9099", "#ffeda0", "#fc4e2a")) +
theme_linedraw() +
theme(axis.text = element_text(color = "black", size = 10),
strip.background = element_blank(),
strip.text.y = element_text(color = "black", size = 11, hjust = 0, vjust = 1),
legend.title = element_blank(),
legend.position = "bottom",
panel.grid.minor = element_blank(),
panel.grid.major = element_line(color = "grey")) +
labs(y = bquote("Cells "~mL^-1), x = "") +
guides(fill=guide_legend(override.aes=list(shape=22)))
# dev.off()Generate complete table with all cell count data
# head(prok)
proks_allyrs_OUTPUT <- prok %>%
separate(SAMPLE, c("VENTSITE", "NAME"), sep = "-", remove = FALSE) %>%
mutate(YEAR = 2020,
Origin = "") %>%
mutate(VENTSITE = case_when(
VENTSITE == "VD" ~ "Von Damm",
TRUE ~ VENTSITE
)) %>%
separate(BAC, c("ORIGSAMPLE", "ID_number"), sep = " ") %>%
select(YEAR, VENTSITE, NAME, REP, CELLML, ORIGSAMPLE, ID_number, Origin) %>%
bind_rows(prok_prev_formatted) %>%
data.frame## Warning: Expected 2 pieces. Additional pieces discarded in 1 rows [27].## Warning: Expected 2 pieces. Missing pieces filled with `NA` in 12 rows [17, 18,
## 19, 20, 21, 33, 34, 35, 36, 39, 40, 41].# View(proks_allyrs_OUTPUT)
# write_delim(proks_allyrs_OUTPUT, path = "input-data/cell-counts-MCR-multiyr.txt", delim = "\t")# head(proks_allyrs_OUTPUT)
# unique(proks_allyrs_OUTPUT$NAME)
# View(proks_allyrs_OUTPUT %>%
# filter(!is.na(as.numeric(CELLML))) %>%
# mutate(TYPE = case_when(
# NAME == "BSW" ~ "Background",
# NAME == "Plume" ~ "Plume",
# TRUE ~ "Vent"
# )) %>%
# mutate(VENTSITE = case_when(
# VENTSITE == "VD" ~ "Von Damm",
# TRUE ~ VENTSITE
# )) %>%
# group_by(VENTSITE, YEAR, TYPE) %>%
# summarise(MEAN = mean(as.numeric(CELLML))) %>%
# pivot_wider(names_from = YEAR, values_from = MEAN, values_fill = NA))
# ?pivot_wider4 Determine grazing effect
Calculate FLP per eukaryotic cell over time. Goal is to make these calculations and then determine best fit line. Slope of best fit line is the grazing rate. Need to take into account euk cells with FLPs and then the euk cells withOUT FLPs, these will be zeroes to take into account for FLPs/euk averages.
Retain previously generated dataframes:
# save(counts_cellsml_all, counts_cellsml_avg, counts_occur, file = "MCR-cellcount-dfs")
load("MCR-cellcount-dfs", verbose = TRUE)## Loading objects:
## counts_cellsml_all
## counts_cellsml_avg
## counts_occur# Euk cell counts, all and averaged
head(counts_cellsml_all)## SAMPLE Site Name EXPID TimePoint Replicate VOL
## 1 Piccard-LotsOShrimp Piccard LotsOShrimp T0-Rep3 T0 Rep3 50
## 2 Piccard-LotsOShrimp Piccard LotsOShrimp T15-Rep3 T15 Rep3 50
## 3 Piccard-LotsOShrimp Piccard LotsOShrimp T20-Rep3 T20 Rep3 50
## 4 Piccard-LotsOShrimp Piccard LotsOShrimp T40-Rep3 T40 Rep3 40
## 5 Piccard-Plume Piccard Plume T0-Rep1 T0 Rep1 150
## 6 Piccard-Plume Piccard Plume T0-Rep2 T0 Rep2 150
## totalFOV nanoAvg nanoVar nanoSd microAvg microVar microSd
## 1 30 0.3666667 0.3781609 1.2298958 0.00000000 0.00000000 0.0000000
## 2 30 0.3000000 0.2172414 0.9321832 0.06666667 0.06436782 0.5074163
## 3 30 0.3333333 0.2988506 1.0933445 0.00000000 0.00000000 0.0000000
## 4 30 0.1666667 0.1436782 0.7580980 0.00000000 0.00000000 0.0000000
## 5 30 0.3000000 0.2862069 1.0699662 0.03333333 0.03333333 0.3651484
## 6 30 0.2000000 0.1655172 0.8136762 0.06666667 0.06436782 0.5074163
## euksAvg euksVar euksSd nanoCONC microCONC eukCONC
## 1 0.3666667 0.3781609 1.229896 230.90630 0.00000 230.90630
## 2 0.3666667 0.2402299 0.980265 188.92333 41.98296 230.90630
## 3 0.3333333 0.2988506 1.093345 209.91481 0.00000 209.91481
## 4 0.1666667 0.1436782 0.758098 131.19676 0.00000 131.19676
## 5 0.3333333 0.3678161 1.212957 62.97444 6.99716 69.97160
## 6 0.2666667 0.2712644 1.041661 41.98296 13.99432 55.97728head(counts_cellsml_avg)## Site Name TimePoint EXP_TYPE IGT_REP VARIABLE MEAN VAR SD SEM
## 1 Piccard LotsOShrimp T0 Bag Bag eukCONC 230.90630 NA NA NA
## 2 Piccard LotsOShrimp T0 Bag Bag microCONC 0.00000 NA NA NA
## 3 Piccard LotsOShrimp T0 Bag Bag nanoCONC 230.90630 NA NA NA
## 4 Piccard LotsOShrimp T15 Bag Bag eukCONC 230.90630 NA NA NA
## 5 Piccard LotsOShrimp T15 Bag Bag microCONC 41.98296 NA NA NA
## 6 Piccard LotsOShrimp T15 Bag Bag nanoCONC 188.92333 NA NA NA
## MIN MAX
## 1 230.90630 230.90630
## 2 0.00000 0.00000
## 3 230.90630 230.90630
## 4 230.90630 230.90630
## 5 41.98296 41.98296
## 6 188.92333 188.92333head(counts_occur)## DATE SAMPLE EXPID VOL MAG FOV nanoNoFLP microNoFLP nanoFLP microFLP
## 1 4/2/20 VD-Plume T0-Rep1 150 100 1 0 0 0 0
## 2 4/2/20 VD-Plume T0-Rep1 150 100 2 0 0 0 0
## 3 4/2/20 VD-Plume T0-Rep1 150 100 3 3 0 0 0
## 4 4/2/20 VD-Plume T0-Rep1 150 100 4 1 0 4 0
## 5 4/2/20 VD-Plume T0-Rep1 150 100 5 0 1 0 0
## 6 4/2/20 VD-Plume T0-Rep1 150 100 6 0 0 0 0
## NOTES DateCompiled nanoFLP_occur microFLP_occur nanoTOTAL
## 1 0 6/7/20 0 0 0
## 2 0 6/7/20 0 0 0
## 3 0 6/7/20 0 0 3
## 4 0 6/7/20 1 0 2
## 5 elongated cell? 6/7/20 0 0 0
## 6 0 6/7/20 0 0 0
## microTOTAL euksTOTAL
## 1 0 0
## 2 0 0
## 3 0 3
## 4 0 2
## 5 1 1
## 6 0 0Supplementary plot showing the concentration of eukaryotic cells over time for each experiment.
# head(counts_cellsml_avg)
vent_ids <- c("BSW","Plume", "Shrimpocalypse", "LotsOShrimp", "X18", "OMT", "Rav2", "MustardStand", "ShrimpHole")
vent_fullname <- c("Background","Plume", "Shrimpocalypse", "Lots 'O Shrimp", "X-18", "Old Man Tree", "Ravelin #2", "Mustard Stand", "Shrimp Hole")
site_ids <- c("VD", "Piccard")
site_fullname <- c("Von Damm", "Piccard")
counts_cellsml_avg$SiteOrder <- factor(counts_cellsml_avg$Site, levels = site_ids, labels = site_fullname)
counts_cellsml_avg$NameOrder <- factor(counts_cellsml_avg$Name, levels = vent_ids, labels = vent_fullname)
# Plot trend line of euk cell count for all experiments
counts_cellsml_avg %>%
filter(VARIABLE == "eukCONC") %>%
unite("Experiment", NameOrder, IGT_REP, sep = "-", remove = FALSE) %>%
ggplot(aes(x = TimePoint, y = MEAN, shape = EXP_TYPE, fill = NameOrder)) +
geom_path(aes(group = Experiment)) +
# geom_errorbar(aes(ymax = (MEAN + SD), ymin = (MEAN - SD)), width = 0.2) +
geom_errorbar(aes(ymax = (MEAN + SEM), ymin = (MEAN - SEM)), width = 0.2) +
geom_point(stat = "identity", size = 2, aes(shape = EXP_TYPE)) +
scale_shape_manual(values = c(21, 24)) +
scale_fill_brewer(palette = "Paired") +
scale_y_log10() +
facet_wrap(SiteOrder ~ EXP_TYPE, scales = "free") +
theme_classic() + theme(strip.background = element_blank(),
legend.title = element_blank(),
title = element_text(size = 7, face = "bold"),
axis.title = element_text(size = 9)) +
labs(title = "Total euk cell counts for each experiment", y = bquote("Average eukaryote cells "~mL^-1), x = "Time point") +
guides(fill=guide_legend(override.aes=list(shape=21)))
## Determine FLP per euk Isolate the euk cell counts with FLPs, these are comma separated and must be separated into rows. Using counts_occur data frame from above.
# Select nano and micro counts with FLPs
counts_sepflp <- counts_occur %>%
filter(!NOTES == "Discard") %>%
filter(!(NOTES == "DTAF stain prevented counts of FLP, Euks only")) %>%
select(DATE, SAMPLE, EXPID, VOL, MAG, FOV, nanoFLP, microFLP) %>%
# Inputs that are comma separated will be split into a new row
separate_rows(microFLP, sep = ",", convert = TRUE) %>%
separate_rows(nanoFLP, sep = ",", convert = TRUE) %>%
# Replace NAs with zeroes
replace_na(list(microFLP = 0, nanoFLP = 0)) %>%
data.frame
## Check, see FOV 23, separated into rows.
# View(counts_sepflp %>%
# filter(SAMPLE == "VD-Rav2" & EXPID == "T10-Rep1"))
# View(counts_occur %>%
# filter(SAMPLE == "VD-Rav2" & EXPID == "T10-Rep1"))Subset counts that are greater than 0, so only euk cells observed to have FLPs. Create new column that calculates FLP per Euk - by dividing by 1.
counts_flp <- counts_sepflp %>%
select(SAMPLE, EXPID, nano_size = nanoFLP, micro_size = microFLP) %>%
pivot_longer(cols = ends_with("_size"), names_to = "SizeFrac", values_to = "num_of_FLP") %>%
filter(num_of_FLP > 0) %>%
separate(SAMPLE, c("Site", "Name"), sep = "-", remove = FALSE) %>%
separate(EXPID, c("TimePoint", "Replicate"), sep = "-", remove = FALSE) %>%
mutate(EXP_TYPE = case_when(
grepl("IGT", Replicate) ~ "IGT",
grepl("Rep", Replicate) ~ "Bag"
)) %>%
mutate(IGT_REP = case_when(
EXP_TYPE == "IGT" ~ Replicate,
EXP_TYPE == "Bag" ~ "Bag")) %>%
group_by(SAMPLE, EXPID, EXP_TYPE, IGT_REP, SizeFrac) %>%
summarise(total_FLP = sum(num_of_FLP),
total_euks_wflp = n(),
.groups = "rowwise") %>%
data.frameOUTPUT COLUMNS: (1) total_FLP = sum of FLPs found inside a euk cell (2) total_euks_wflp = number of euks counted with ingested FLP
Repeat above operation for euk cells without any FLP. Here, subset total number of observations where there was a euk cell without FLP. These need to be counted as euk cell without an FLP.
Below code repeats process and compiles with other FLP/euk cell data.
counts_flp_compiled <- counts_occur %>%
filter(!(NOTES == "Discard")) %>% #Discard bad counts
filter(!(NOTES == "DTAF stain prevented counts of FLP, Euks only")) %>%
type.convert(as.is = TRUE) %>% #modify str() for columns
select(SAMPLE, EXPID, nano_size = nanoNoFLP, micro_size = microNoFLP) %>% #select non flp
pivot_longer(cols = ends_with("_size"), names_to = "SizeFrac", values_to = "num_of_euks") %>%
separate(SAMPLE, c("Site", "Name"), sep = "-", remove = FALSE) %>%
separate(EXPID, c("TimePoint", "Replicate"), sep = "-", remove = FALSE) %>%
mutate(EXP_TYPE = case_when(
grepl("IGT", Replicate) ~ "IGT",
grepl("Rep", Replicate) ~ "Bag"
)) %>%
mutate(IGT_REP = case_when(
EXP_TYPE == "IGT" ~ Replicate,
EXP_TYPE == "Bag" ~ "Bag")) %>%
# filter(num_of_euks > 0) %>% # Remove observed zero counts
group_by(SAMPLE, EXPID, EXP_TYPE, IGT_REP, SizeFrac) %>%
summarise(total_euks_noFLP = sum(num_of_euks),
.groups = "rowwise") %>%
# Join with FLP count information
## SAMPLE, EXPID, EXPTYPE, IGTREP, and SizeFrac variables should match
left_join(counts_flp) %>% # Join with the counts of FLP per euk cell
replace_na(list(total_FLP = 0, total_euks_wflp = 0)) %>% #Replace NAs with zero
data.frame## Joining, by = c("SAMPLE", "EXPID", "EXP_TYPE", "IGT_REP", "SizeFrac")## Check data frames
# View(counts_flp_compiled)
# View(counts_occur)
# (filter(counts_flp_compiled, SAMPLE == "VD-X18", EXPID == "T40-Rep2"))
# (filter(counts_flp_compiled, SAMPLE == "Piccard-Plume", EXPID == "T10-Rep1"))
# head(counts_flp_compiled)Extract total euk value by adding across nano and micro, then combine back with nano and micro counts.
counts_flp_compiled_all <- counts_flp_compiled %>%
# Exclude size fraction:
group_by(SAMPLE, EXPID, EXP_TYPE, IGT_REP) %>%
summarise(total_euks_noFLP = sum(total_euks_noFLP),
total_FLP = sum(total_FLP),
total_euks_wflp = sum(total_euks_wflp),
.groups = "rowwise") %>%
add_column(SizeFrac = "total_euks") %>% #Add SizeFrac column
bind_rows(counts_flp_compiled) %>% # Combine back with flp compiled list
data.frame
# head(counts_flp_compiled_all)
# filter(counts_flp_compiled_all, SAMPLE == "VD-X18", EXPID == "T40-Rep2")
# filter(counts_flp_compiled_all, SAMPLE == "Piccard-Plume", EXPID == "T10-Rep1")All three size fractions are represented, total, micro, and nano.
One outstanding question will be that there were many more nano size grazers than micro sized. Here I am preferring to use the total value (micro + nano), but will need to defend that this is ok.
# unique(counts_flp_compiled_all$SizeFrac)
# unique(counts_flp_compiled_all$SAMPLE)
# head(counts_flp_compiled_all[1:2,])4.1 Calculate FLP per euk cell calculation
Import metadata which has the extact minutes for each time point used.
metadata <- read.delim("input-data/flp-exp-metadata-compiled.txt")
head(metadata)## SAMPLE TimePoint Minutes SiteOrigin REP
## 1 VD-MustardStand T0 0 Vent Bag
## 2 VD-MustardStand T10 10 Vent Bag
## 3 VD-MustardStand T15 15 Vent Bag
## 4 VD-MustardStand T20 20 Vent Bag
## 5 VD-MustardStand T40 40 Vent Bag
## 6 Piccard-Plume T0 0 Plume BagAdd metadata information to FLP data, reformat sample names, and calculate FLP per euk as the total FLP divided by the total number of euk cells counted.
counts_flp_calcs_all <- counts_flp_compiled_all %>%
# Add in metadata
# IGTXb are replicate counts, include them as replicates!
separate(EXPID, c("TimePoint", "REP"), sep = "-", remove = FALSE) %>% mutate(
REP = ifelse(grepl("IGT5b", REP), "IGT5", REP),
REP = ifelse(grepl("IGT4b", REP), "IGT4", REP),
REP = ifelse(grepl("Bag", EXP_TYPE), "Bag", REP)) %>%
left_join(metadata, by = c("SAMPLE" = "SAMPLE", "TimePoint" = "TimePoint", "REP" = "REP")) %>%
separate(SAMPLE, c("Site", "Name"), sep = "-", remove = FALSE) %>%
separate(EXPID, c("TimePoint", "Replicate_ID"), sep = "-", remove = FALSE) %>%
## Treat repeated IGT counts completely separate
# group_by(SAMPLE, Site, Name, EXPID, TimePoint, Replicate_ID, EXP_TYPE, IGT_REP, SizeFrac) %>%
## Treat repeated IGT counts as replicates (e.g., IGT4b and IGT4 == IGT4)
group_by(SAMPLE, Site, Name, EXPID, TimePoint, Replicate_ID, EXP_TYPE, REP, SizeFrac) %>%
# FLPperEuk is the total FLP divided by the total number of euk cells counted
mutate(FLPperEuk = total_FLP/(sum(total_euks_noFLP, total_euks_wflp))) %>%
unite("Experiment", Name, REP, sep = "-", remove = FALSE) %>%
data.frame
# filter(counts_flp_calcs_all, SAMPLE == "VD-X18", EXPID == "T40-Rep2")
# filter(counts_flp_calcs_all, SAMPLE == "Piccard-Plume", EXPID == "T10-Rep1")View output
# head(counts_flp_calcs_all[1:2,])
unique(counts_flp_calcs_all$Experiment)## [1] "LotsOShrimp-Bag" "Plume-Bag" "Shrimpocalypse-IGT3"
## [4] "Shrimpocalypse-Bag" "BSW-Bag" "MustardStand-Bag"
## [7] "OMT-IGT4" "Rav2-IGT4" "Rav2-IGT5"
## [10] "Rav2-Bag" "ShrimpHole-Bag" "X18-Bag"# View(counts_flp_calcs_all)COLS: Timepoint, Minutes = time point label, actual incubated minutes
COLS: Replicate_ID, REP, and IGT_REP = full replicate identified for IGTs and Bags, designation of biological replicates, and designation of technical replicates for IGT experiments
4.2 Calculate linear regression to obtain slope
Use lm() function in R to calculate linear regression for each experiment. Slope equates to grazing rate. Function inputs the FLP per euk cell data, performs regression and then adds a column for slope and r-squared values.
calculate_lm <- function(df){
regression_1 <- df %>%
type.convert(as.is = TRUE) %>%
## Keep technical replicates separate for IGTs
# group_by(SAMPLE, Site, Experiment, Name, IGT_REP, SizeFrac) %>%
# nest(-SAMPLE, -Site, -Experiment, -Name, -IGT_REP, -SizeFrac) %>%
## Combine technical replicates for IGTs
group_by(SAMPLE, Site, Experiment, Name, REP, SizeFrac) %>%
nest(-SAMPLE, -Site, -Experiment, -Name, -REP, -SizeFrac) %>%
mutate(lm_fit = map(data, ~lm(FLPperEuk ~ Minutes, data = .)),
tidied = map(lm_fit, tidy)) %>%
unnest(tidied) %>%
# select(SAMPLE, Site, Experiment, Name, IGT_REP, SizeFrac, term, estimate) %>%
select(SAMPLE, Site, Experiment, Name, REP, SizeFrac, term, estimate) %>%
pivot_wider(names_from = term, values_from = estimate) %>%
data.frame
# colnames(regression_1) <- c("SAMPLE", "Site", "Experiment", "Name", "IGT_REP",
# "SizeFrac", "INTERCEPT", "SLOPE")
colnames(regression_1) <- c("SAMPLE", "Site",
"Experiment", "Name", "REP",
"SizeFrac", "INTERCEPT", "SLOPE")
out_regression <- df %>%
# group_by(SAMPLE, Site, Experiment, Name, IGT_REP, SizeFrac) %>%
# nest(-SAMPLE, -Site, -Experiment, -Name, -IGT_REP, -SizeFrac) %>%
group_by(SAMPLE, Site, Experiment, Name, REP, SizeFrac) %>%
nest(-SAMPLE, -Site, -Experiment, -Name, -REP, -SizeFrac) %>%
mutate(lm_fit = map(data, ~lm(FLPperEuk ~ Minutes, data = .)),
glanced = map(lm_fit, glance)) %>%
unnest(glanced) %>%
# select(SAMPLE, Site, Experiment, Name, IGT_REP, SizeFrac, r.squared) %>%
select(SAMPLE, Site, Experiment, Name, REP, SizeFrac, r.squared) %>%
right_join(regression_1) %>%
right_join(df) %>%
data.frame
out_regression$SITE <- factor(out_regression$Site, levels = c("VD", "Piccard"))
out_regression$TYPE <- factor(out_regression$EXP_TYPE, levels = c("Bag", "IGT"))
return(out_regression)
}Note that an error may occur when running the below function. This is due to the fact that some experiments did not have replicates.
calcs_wslope_regression <- calculate_lm(counts_flp_calcs_all)## Warning: All elements of `...` must be named.
## Did you want `data = c(EXPID, TimePoint, Replicate_ID, EXP_TYPE, IGT_REP, total_euks_noFLP,
## total_FLP, total_euks_wflp, Minutes, SiteOrigin, FLPperEuk)`?## Warning in summary.lm(x): essentially perfect fit: summary may be unreliable## Warning: All elements of `...` must be named.
## Did you want `data = c(EXPID, TimePoint, Replicate_ID, EXP_TYPE, IGT_REP, total_euks_noFLP,
## total_FLP, total_euks_wflp, Minutes, SiteOrigin, FLPperEuk)`?## Warning in summary.lm(x): essentially perfect fit: summary may be unreliable
## Warning in summary.lm(x): essentially perfect fit: summary may be unreliable## Joining, by = c("SAMPLE", "Site", "Experiment", "Name", "REP", "SizeFrac")
## Joining, by = c("SAMPLE", "Site", "Experiment", "Name", "REP", "SizeFrac")# head(calcs_wslope_regression[1:2,])
# View(calcs_wslope_regression)Use below commented out to recalculate one linear regression. Above function uses the nest() capability of tidyverse. Below, one experiment is subset to check the value.
# Extract only plume-bag experiment from VD
# tmp_plume <- filter(counts_flp_calcs_all, Experiment == "Plume-Bag") %>% filter(Site == "VD") %>% filter(SizeFrac == "total_euks")
# tmp_plume # View
# Perform linear regression
# lm_out <- lm(FLPperEuk ~ Minutes, data = tmp_plume)
# Check output
# summary(lm_out)
# lm_out$coefficients #Intercept=intercept #Minutes = SLOPE
# Compare with nested function output
# filter(calcs_wslope_regression, Experiment == "Plume-Bag") %>% filter(Site == "VD") %>% filter(SizeFrac == "total_euks") %>% head# View(calcs_wslope_regression)
# head(calcs_wslope_regression)
# unique(calcs_wslope_regression$Site)4.3 Plot linear regression trend
calcs_wslope_regression %>%
filter(SizeFrac == "total_euks") %>%
# filter(TYPE != "IGT") %>%
unite(EXPERIMENT, SITE, Experiment, sep = " ", remove = FALSE) %>%
ggplot(aes(x = Minutes, y = FLPperEuk, fill = Site, shape = TYPE, group = Experiment)) +
geom_abline(aes(slope = SLOPE, intercept = INTERCEPT), color = "black", linetype = "dashed", size = 1) +
geom_point(stat = "identity", color = "black",
size = 2, aes(shape = TYPE, fill = Site)) +
scale_shape_manual(values = c(21, 24)) +
scale_fill_manual(values = c("#de2d26", "#1c9099")) +
labs(x = "Minutes", y = bquote("FLP"~eukaryote^-1), title = "Grazing experiment regression") +
facet_wrap(. ~ EXPERIMENT) +
# Report r.squared
geom_text(aes(x = 42, y = max(FLPperEuk), label = paste(round(SLOPE, 4))),
vjust = 1, hjust = 0, size = 3) +
theme_bw() +
theme(strip.background = element_blank(),
strip.text = element_text(color = "black", size = 7),
legend.title = element_blank(),
legend.position = "right")
Data points represent the FLP per euk cells (based on total eukaryote cells counts). Y-axis represents the duration of incubation (in minutes). The dashed purple line reprents the slope and intercept of the experiment.
4.4 Remove Tf-IGT experiments
IGT experiment results appear to have bottle effect, especially in the final time point. Additionally, due to the lack of biological replicates in the IGT experiments, technical replicates are treated as biological replicates in the regression below.
IGT_lm_woTf <- counts_flp_calcs_all %>%
# Select only IGT experiments with total eukaryotes, remove Tf (T3)
filter(SizeFrac == "total_euks") %>%
filter(EXP_TYPE == "IGT" & !(TimePoint == "T3")) %>%
add_column(IGT_cor = "rm Tf") %>%
data.frame
# Recalculate lm(), keep replicates separate
igt_regression_noTf <- calculate_lm(IGT_lm_woTf) # Recalculate## Warning: All elements of `...` must be named.
## Did you want `data = c(EXPID, TimePoint, Replicate_ID, EXP_TYPE, IGT_REP, total_euks_noFLP,
## total_FLP, total_euks_wflp, Minutes, SiteOrigin, FLPperEuk,
## IGT_cor)`?
## Warning: All elements of `...` must be named.
## Did you want `data = c(EXPID, TimePoint, Replicate_ID, EXP_TYPE, IGT_REP, total_euks_noFLP,
## total_FLP, total_euks_wflp, Minutes, SiteOrigin, FLPperEuk,
## IGT_cor)`?## Joining, by = c("SAMPLE", "Site", "Experiment", "Name", "REP", "SizeFrac")
## Joining, by = c("SAMPLE", "Site", "Experiment", "Name", "REP", "SizeFrac")Plot IGT grazing experiment with newly calculated grazing effect.
igt_regression_noTf %>%
# filter(SizeFrac == "total_euks") %>%
# filter(TYPE != "IGT") %>%
unite(EXPERIMENT, SITE, Experiment, sep = " ", remove = FALSE) %>%
ggplot(aes(x = Minutes, y = FLPperEuk, fill = Site, shape = TYPE, group = Experiment)) +
geom_abline(aes(slope = SLOPE, intercept = INTERCEPT), color = "black", linetype = "dashed", size = 1) +
geom_point(stat = "identity", color = "black",
size = 2, aes(shape = TYPE, fill = Site)) +
scale_shape_manual(values = c(24)) +
scale_fill_manual(values = c("#de2d26", "#1c9099")) +
labs(x = "Minutes", y = bquote("FLP"~eukaryote^-1), title = "Grazing experiment regression") +
facet_wrap(. ~ EXPERIMENT) +
# Report r.squared
geom_text(aes(x = 5, y = max(FLPperEuk), label = paste(round(SLOPE, 4))),
vjust = 1, hjust = 0, size = 3) +
theme_bw() +
theme(strip.background = element_blank(),
strip.text = element_text(color = "black", size = 7),
legend.title = element_blank(),
legend.position = "right")
IGT grazing regression without the Tf data point appear more consistent. Keeping that one.
4.5 Compile grazing experiment results
calcs_wslope_regression_update <- calcs_wslope_regression %>%
filter(TYPE != "IGT") %>%
bind_rows(igt_regression_noTf %>% select(-IGT_cor)) %>%
data.frame
# Factor
vent_ids <- c("BSW","Plume", "Shrimpocalypse", "LotsOShrimp", "X18", "OMT", "Rav2", "MustardStand", "ShrimpHole")
vent_fullname <- c("Background","Plume", "Shrimpocalypse", "Lots 'O Shrimp", "X-18", "Old Man Tree", "Ravelin #2", "Mustard Stand", "Shrimp Hole")
site_ids <- c("VD", "Piccard")
site_fullname <- c("Von Damm", "Piccard")
# Factor for shipboard
calcs_wslope_regression_update$SiteOrder <- factor(calcs_wslope_regression_update$Site, levels = site_ids, labels = site_fullname)
calcs_wslope_regression_update$NameOrder <- factor(calcs_wslope_regression_update$Name, levels = vent_ids, labels = vent_fullname)
# dim(calcs_wslope_regression); dim(calcs_wslope_regression_update)
head(calcs_wslope_regression_update)## SAMPLE Site Experiment Name REP SizeFrac
## 1 Piccard-LotsOShrimp Piccard LotsOShrimp-Bag LotsOShrimp Bag total_euks
## 2 Piccard-LotsOShrimp Piccard LotsOShrimp-Bag LotsOShrimp Bag total_euks
## 3 Piccard-LotsOShrimp Piccard LotsOShrimp-Bag LotsOShrimp Bag total_euks
## 4 Piccard-LotsOShrimp Piccard LotsOShrimp-Bag LotsOShrimp Bag total_euks
## 5 Piccard-Plume Piccard Plume-Bag Plume Bag total_euks
## 6 Piccard-Plume Piccard Plume-Bag Plume Bag total_euks
## r.squared INTERCEPT SLOPE EXPID TimePoint Replicate_ID EXP_TYPE
## 1 0.66025258 0.3517002 -0.007605829 T0-Rep3 T0 Rep3 Bag
## 2 0.66025258 0.3517002 -0.007605829 T15-Rep3 T15 Rep3 Bag
## 3 0.66025258 0.3517002 -0.007605829 T20-Rep3 T20 Rep3 Bag
## 4 0.66025258 0.3517002 -0.007605829 T40-Rep3 T40 Rep3 Bag
## 5 0.02810998 0.7435088 0.005362856 T0-Rep1 T0 Rep1 Bag
## 6 0.02810998 0.7435088 0.005362856 T0-Rep2 T0 Rep2 Bag
## IGT_REP total_euks_noFLP total_FLP total_euks_wflp Minutes SiteOrigin
## 1 Bag 9 3 2 0 Vent
## 2 Bag 8 4 3 15 Vent
## 3 Bag 9 2 1 20 Vent
## 4 Bag 5 0 0 40 Vent
## 5 Bag 6 6 4 0 Plume
## 6 Bag 5 5 3 0 Plume
## FLPperEuk SITE TYPE SiteOrder NameOrder
## 1 0.2727273 Piccard Bag Piccard Lots 'O Shrimp
## 2 0.3636364 Piccard Bag Piccard Lots 'O Shrimp
## 3 0.2000000 Piccard Bag Piccard Lots 'O Shrimp
## 4 0.0000000 Piccard Bag Piccard Lots 'O Shrimp
## 5 0.6000000 Piccard Bag Piccard Plume
## 6 0.6250000 Piccard Bag Piccard Plume4.6 Check FLP control experiments
All incubations had control experiments run alongside them. This was to ensure added FLP did not decrease or change in concentration over time.
bac_ctrl <- read.delim("input-data/bac-counts-compiled.txt")
# dim(bac_ctrl)
dtaf <- bac_ctrl %>%
separate(SampleID, c("exp", "Replicate", "TimePoint"), sep = "-", remove = FALSE) %>%
separate(Site, c("Site", "Name"), sep = "-", remove = FALSE) %>%
filter(Stain == "DTAF") %>%
data.frame## Warning: Expected 2 pieces. Additional pieces discarded in 17 rows [33, 34, 35,
## 36, 37, 38, 39, 40, 41, 42, 43, 44, 45, 46, 47, 48, 49].# View(bac_ctrl)
# head(dtaf)
dtaf_avg <- dtaf %>%
group_by(TimePoint, Stain, Site, Name) %>%
summarise(Avg_cellsperml = mean(Cells.ml)) %>%
data.frame## `summarise()` has grouped output by 'TimePoint', 'Stain', 'Site'. You can override using the `.groups` argument.dtaf_avg %>%
filter(Site != "IGT") %>%
ggplot(aes(x = TimePoint, y = Avg_cellsperml, fill = Name, shape = Site)) +
geom_rect(data = filter(dtaf_avg, TimePoint == "T0", Site != "IGT"), aes(
ymin = (Avg_cellsperml-(0.1*Avg_cellsperml)),
ymax = (Avg_cellsperml+(0.1*Avg_cellsperml))), color = NA, alpha = 0.4, xmin = 0, xmax = 6, fill = "black") +
geom_line(aes(group = Name)) +
geom_point(stat = "identity", aes(shape = Site, fill = Name), size = 2) +
# scale_fill_manual(values = c("black","#9970ab", "#5aae61")) +
facet_wrap(Name ~ Site) +
scale_y_log10() +
theme_bw() + theme(strip.background = element_blank(),
legend.title = element_blank(),
axis.text = element_text(size = 10, color = "black"),
title = element_text(size = 10, face = "bold"),
axis.title = element_text(size = 9)) +
labs(title = "FLP counts for controls", y = bquote("Log FLP "~mL^-1), x = "Time point")
Control experiments look ok for now.
dtaf_avg %>%
filter(Site == "IGT") %>%
ggplot(aes(x = TimePoint, y = Avg_cellsperml, fill = Name, shape = Site)) +
geom_rect(data = filter(dtaf_avg, TimePoint == "T0", Site == "IGT"), aes(
ymin = (Avg_cellsperml-(0.1*Avg_cellsperml)),
ymax = (Avg_cellsperml+(0.1*Avg_cellsperml))), color = NA, alpha = 0.4, xmin = 0, xmax = 6, fill = "black") +
geom_line(aes(group = Name)) +
geom_point(stat = "identity", aes(shape = Site, fill = Name), size = 2) +
# scale_fill_manual(values = c("black","#9970ab", "#5aae61")) +
facet_wrap(Name ~ Site) +
scale_y_log10() +
theme_bw() + theme(strip.background = element_blank(),
legend.title = element_blank(),
axis.text = element_text(size = 10, color = "black"),
title = element_text(size = 10, face = "bold"),
axis.title = element_text(size = 9)) +
labs(title = "FLP counts for controls", y = bquote("Log FLP "~mL^-1), x = "Time point")
4.7 Table with grazing values
head(calcs_wslope_regression_update)## SAMPLE Site Experiment Name REP SizeFrac
## 1 Piccard-LotsOShrimp Piccard LotsOShrimp-Bag LotsOShrimp Bag total_euks
## 2 Piccard-LotsOShrimp Piccard LotsOShrimp-Bag LotsOShrimp Bag total_euks
## 3 Piccard-LotsOShrimp Piccard LotsOShrimp-Bag LotsOShrimp Bag total_euks
## 4 Piccard-LotsOShrimp Piccard LotsOShrimp-Bag LotsOShrimp Bag total_euks
## 5 Piccard-Plume Piccard Plume-Bag Plume Bag total_euks
## 6 Piccard-Plume Piccard Plume-Bag Plume Bag total_euks
## r.squared INTERCEPT SLOPE EXPID TimePoint Replicate_ID EXP_TYPE
## 1 0.66025258 0.3517002 -0.007605829 T0-Rep3 T0 Rep3 Bag
## 2 0.66025258 0.3517002 -0.007605829 T15-Rep3 T15 Rep3 Bag
## 3 0.66025258 0.3517002 -0.007605829 T20-Rep3 T20 Rep3 Bag
## 4 0.66025258 0.3517002 -0.007605829 T40-Rep3 T40 Rep3 Bag
## 5 0.02810998 0.7435088 0.005362856 T0-Rep1 T0 Rep1 Bag
## 6 0.02810998 0.7435088 0.005362856 T0-Rep2 T0 Rep2 Bag
## IGT_REP total_euks_noFLP total_FLP total_euks_wflp Minutes SiteOrigin
## 1 Bag 9 3 2 0 Vent
## 2 Bag 8 4 3 15 Vent
## 3 Bag 9 2 1 20 Vent
## 4 Bag 5 0 0 40 Vent
## 5 Bag 6 6 4 0 Plume
## 6 Bag 5 5 3 0 Plume
## FLPperEuk SITE TYPE SiteOrder NameOrder
## 1 0.2727273 Piccard Bag Piccard Lots 'O Shrimp
## 2 0.3636364 Piccard Bag Piccard Lots 'O Shrimp
## 3 0.2000000 Piccard Bag Piccard Lots 'O Shrimp
## 4 0.0000000 Piccard Bag Piccard Lots 'O Shrimp
## 5 0.6000000 Piccard Bag Piccard Plume
## 6 0.6250000 Piccard Bag Piccard Plume# Generate final table
table_grazerate <- calcs_wslope_regression_update %>%
filter(SizeFrac == "total_euks") %>%
group_by(SAMPLE, SiteOrder, NameOrder, SLOPE, EXP_TYPE) %>%
mutate(No_of_Replicates = n_distinct(Replicate_ID)) %>%
select(SAMPLE, Site=SiteOrder, Name=NameOrder, RATE = SLOPE, EXP_TYPE, No_of_Replicates) %>%
distinct() %>%
data.frame
# View(table_grazerate)
# write_delim(table_grazerate, path = "Table-grazingrates.txt", delim = "\t")library(gt)
table_grazerate %>%
group_by(Site, EXP_TYPE) %>%
gt() #%>% | SAMPLE | Name | RATE | No_of_Replicates |
|---|---|---|---|
| Piccard - Bag | |||
| Piccard-LotsOShrimp | Lots 'O Shrimp | -0.007605829 | 1 |
| Piccard-Plume | Plume | 0.005362856 | 3 |
| Piccard-Shrimpocalypse | Shrimpocalypse | 0.015686872 | 1 |
| Von Damm - Bag | |||
| VD-BSW | Background | 0.002958889 | 3 |
| VD-MustardStand | Mustard Stand | -0.005445545 | 2 |
| VD-Plume | Plume | 0.005274231 | 3 |
| VD-Rav2 | Ravelin #2 | 0.003470217 | 3 |
| VD-ShrimpHole | Shrimp Hole | -0.001967253 | 2 |
| VD-X18 | X-18 | 0.001744429 | 2 |
| Piccard - IGT | |||
| Piccard-Shrimpocalypse | Shrimpocalypse | 0.016794258 | 1 |
| Von Damm - IGT | |||
| VD-OMT | Old Man Tree | 0.014510943 | 2 |
| VD-Rav2 | Ravelin #2 | 0.015395240 | 2 |
| VD-Rav2 | Ravelin #2 | 0.001603035 | 2 |
# fmt_number(columns = 3, decimals = 4)Add FLP conc to table
# head(table_grazerate)
dtaf_igt <- 5352.8278 # Manually insert FLP concentration for IGT experiments; this value is estimated from how IGT FLP spike-ins were calculated
table_grazerate_wflp <- bac_ctrl %>%
filter(FLP_t0 == "use") %>%
add_column(EXP_TYPE = "Bag") %>%
group_by(Site, EXP_TYPE) %>%
summarise(FLP_conc = mean(Cells.ml)) %>%
right_join(table_grazerate, by = c("Site" = "SAMPLE", "EXP_TYPE" = "EXP_TYPE")) %>%
mutate(FLP_conc = ifelse(EXP_TYPE == "IGT", dtaf_igt, FLP_conc)) %>%
data.frame## `summarise()` has grouped output by 'Site'. You can override using the `.groups` argument.# View(table_grazerate_wflp)Negative grazing rates equal undetected (change to zero).
# head(table_grazerate)
bsw <- c("Plume", "Background")
grazerate_plot <- table_grazerate %>%
type.convert(as.is = TRUE) %>%
mutate(detected = case_when(
RATE < 0 ~ "Not detected",
TRUE ~ "Detected")) %>%
mutate(type = case_when(
Name %in% bsw ~ Name,
TRUE ~ paste("Vent", EXP_TYPE, sep="-")
)) %>%
mutate(GRAZE_RATE = case_when(
RATE < 0 ~ 0,
TRUE ~ RATE
)) %>%
mutate(type_site = case_when(
Name %in% bsw ~ Name,
TRUE ~ "Vent"
)) %>%
# filter(detected == "Detected") %>%
data.frame# Factoring
type_order <- c("Vent-Bag", "Vent-IGT", "Plume", "Background")
grazerate_plot$TYPE <- factor(grazerate_plot$type, levels = type_order)
vent_ids <- c("BSW","Plume", "Shrimpocalypse", "LotsOShrimp", "X18", "OMT", "Rav2", "MustardStand", "ShrimpHole")
vent_fullname <- c("Background","Plume", "Shrimpocalypse", "Lots 'O Shrimp", "X-18", "Old Man Tree", "Ravelin #2", "Mustard Stand", "Shrimp Hole")
site_ids <- c("VD", "Piccard")
site_fullname <- c("Von Damm", "Piccard")
vent_colors <- c("#0868ac", "#41ab5d", "#e7298a", "#c994c7", "#fc4e2a", "#fed976", "#6a51a3", "#ffeda0", "#a1d99b")
names(vent_colors) <- vent_fullname
grazerate_plot$NAME <- factor(grazerate_plot$Name, levels = vent_fullname)# svg("", h =, w = )
grazing_min_plot <- grazerate_plot %>%
ggplot(aes(y = GRAZE_RATE, x = NAME, shape = EXP_TYPE, fill = Site)) +
geom_jitter(stat = "identity", aes(shape = EXP_TYPE, fill = Site),
color = "black", size = 3, width = 0.3) +
scale_shape_manual(values = c(21, 24)) +
scale_fill_manual(values = c("#de2d26", "#1c9099")) +
facet_grid(.~Site, space = "free", scales = "free") +
# coord_flip() +
theme_minimal() +
theme(panel.grid.major = element_line(), panel.grid.minor = element_blank(),
panel.background = element_blank(),
axis.line = element_line(colour = "black"),
axis.text.x = element_text(color="black", size = 12,
angle = 45, hjust = 1, vjust = 1),
axis.text.y = element_text(color="black", size = 12),
axis.title =element_text(color="black", size = 12),
axis.ticks = element_line(),
strip.text =element_blank(), legend.title = element_blank())+
guides(fill = guide_legend(override.aes = list(shape = c(21))),
shape = guide_legend(override.aes = list(fill = "black"))) +
labs(x = "", y = bquote("FLP " ~grazer^-1 ~min^-1))
# dev.off()
grazing_min_plot
5 Plot cell conc & grazing
Generate final plots that show eukaryote and prokaryote biomass and cell concentration. Save final table reporting carbon estimates.
# Subset the average in situ prok cells/ml for non-background samples
tmp <- filter(insitu_proks, Name != "BSW", Name != "Plume") %>% select(MEAN)
avg_insitu <- mean(tmp$MEAN)
# head(insitu_proks)
# Add to master table with data
table_grazerate_wflp_wprok <- insitu_proks %>%
select(Site = SAMPLE, Prok_conc = MEAN, Prok_sem = SEM) %>%
right_join(table_grazerate_wflp) %>%
mutate(Prok_conc = ifelse(is.na(Prok_conc), avg_insitu, Prok_conc)) %>%
data.frame## Joining, by = "Site"Table with all cell count data
# head(table_grazerate_wflp_wprok)
table_grazerate_wflp_wprok_weuk <- plot_euk_format %>%
select(Name = NameOrder, Site.y.y = SiteOrder, euk_conc = avg_conc, EXP_TYPE, euk_conc_sem = SEM_conc) %>%
right_join(table_grazerate_wflp_wprok) %>%
select(SITE = Site.y.y, NAME = Name, EXP = EXP_TYPE, SAMPLE = Site, RATE_min = RATE, REPs = No_of_Replicates, FLP_ml = FLP_conc, PROK_ml = Prok_conc, PROK_sem = Prok_sem, EUK_ml = euk_conc, EUK_sem = euk_conc_sem) %>%
mutate(RATE_min = case_when(
RATE_min < 0 ~ 0,
TRUE ~ RATE_min
)) %>%
data.frame## Joining, by = c("Name", "Site.y.y", "EXP_TYPE")6 Grazing rate calculations
# head(table_grazerate_wflp_wprok_weuk)Based on Unrein et al. 2007, we use the estimated grazing rate, in situ prok abundance, in situ euk abundance, and the concentration of FLP to make additional estimates.
SKH: add citations for the below calulcations.
table_wcalcs <- table_grazerate_wflp_wprok_weuk %>%
# Ingestion rate per hour
mutate(RATE_hr = (RATE_min * 60),
RATE_day = (RATE_hr * 24), #Compare to GR?
# FLP concentration per L
FLP_L = (FLP_ml * 1000),
# mL per grazer per hr
CLEARANCE_RATE_ml = (RATE_hr/FLP_ml),
# nL per grazer per hour
CLEARANCE_RATE_nL = ((RATE_hr/FLP_ml)/1.00E+6),
# proks per grazer per hr
SPEC_GRAZE_RATE_hr = (CLEARANCE_RATE_ml * PROK_ml),
# proks per grazer per day
GRAZE_RATE_DAY = (24 * SPEC_GRAZE_RATE_hr),
# proks per ml per hr
GRAZING_EFFECT_hr = (SPEC_GRAZE_RATE_hr * EUK_ml),
GRAZING_EFFECT_hr_min = (SPEC_GRAZE_RATE_hr * (EUK_ml - EUK_sem)),
GRAZING_EFFECT_hr_max = (SPEC_GRAZE_RATE_hr * (EUK_ml + EUK_sem)),
# cells per ml per day
GRAZING_EFFECT_day = ((SPEC_GRAZE_RATE_hr * 24) * EUK_ml),
# Percentage per day
BAC_TURNOVER_PERC = 100*(GRAZING_EFFECT_day / PROK_ml),
BAC_TURNOVER_PERC_min = 100*(GRAZING_EFFECT_day / (PROK_ml - PROK_sem)),
BAC_TURNOVER_PERC_max = 100*(GRAZING_EFFECT_day / (PROK_ml + PROK_sem))) %>%
data.frame
# View(table_wcalcs)Explanation of units for table with calculated values.
RATE_min & RATE_hr = Grazing rate as ‘FLPs per grazer per minute’ and per hour
CLEARANCE_RATE = ml or nL per grazer per hour
SPEC_GRAZE_RATE (Specific grazing rate) = Prokaryotes per grazer per hour
GRAZING EFFECT = bacteria per ml per hour
Bacterial turnover rate = % per day
7 Carbon biomass table construction
Unrein, F., Massana, R., Alonso-Sáez, L., and Gasol, J.M. (2007) Significant year-round effect of small mixotrophic flagellates on bacterioplankton in an oligotrophic coastal system. Limnol Oceanogr 52: 456–469.
# colnames(table_wcalcs)
bkgd <- c("Background", "Plume")
table_wcalcs %>%
mutate(loc_type = case_when(
NAME %in% bkgd ~ "Background",
TRUE ~ "Vent fluid"
)) %>%
# group_by(loc_type, SITE, EXP) %>%
select(-SAMPLE) %>%
gt(
groupname_col = c("SITE", "EXP", "loc_type"),
rowname_col = "NAME"
) %>%
cols_label(RATE_min = html("minute<sup>-1</sup>"),
RATE_hr = html("hour<sup>-1</sup>"),
RATE_day = html("day<sup>-1</sup>"),
REPs = html("# of incubations"),
FLP_ml = html("FLP ml<sup>-1</sup>"),
PROK_ml = html("Prokaryote cells ml<sup>-1</sup>"),
PROK_sem = html("SEM prokaryote cells ml<sup>-1</sup>"),
EUK_ml = html("Eukaryote cells ml<sup>-1</sup>"),
EUK_sem = html("SEM eukaryote cells ml<sup>-1</sup>"),
FLP_L = html("FLP L<sup>-1</sup>"),
CLEARANCE_RATE_ml = html("ml grazer<sup>-1</sup> hr<sup>-1</sup>"),
CLEARANCE_RATE_nL = html("nl grazer<sup>-1</sup> hr<sup>-1</sup>"),
SPEC_GRAZE_RATE_hr = html("Prokaryote grazer<sup>-1</sup> hr<sup>-1</sup>"),
GRAZE_RATE_DAY = html("Prokaryote grazer<sup>-1</sup> day<sup>-1</sup>"),
GRAZING_EFFECT_hr = html("Prokaryote ml<sup>-1</sup> hr<sup>-1</sup>"),
GRAZING_EFFECT_hr_min = html("MIN"),
GRAZING_EFFECT_hr_max = html("MAX"),
GRAZING_EFFECT_day = html("Prokaryote ml<sup>-1</sup> day<sup>-1</sup>"),
BAC_TURNOVER_PERC = html("Bacteria turnover % day<sup>-1</sup>"),
BAC_TURNOVER_PERC_min = html("MIN"),
BAC_TURNOVER_PERC_max = html("MAX")) %>%
tab_spanner(
label = (html("Turnover")),
columns = starts_with("BAC_TURNOVER")
) %>%
tab_spanner(
label = (html("Grazing rate: prokaryote cells consumed")),
columns = starts_with("GRAZING_EFFECT")
) %>%
tab_spanner(
label = (html("ml grazer<sup>-1</sup> hr<sup>-1</sup>")),
columns = c(CLEARANCE_RATE_ml, CLEARANCE_RATE_nL)
) %>%
tab_spanner(
label = html("Specific grazing rate"),
columns = c(SPEC_GRAZE_RATE_hr, GRAZE_RATE_DAY)
) %>%
tab_spanner(
label = (html("FLPs grazer<sup>-1</sup>")),
columns = c(RATE_hr, RATE_min, RATE_day),
) %>%
tab_spanner(
label = (html("Cell counts")),
columns = c(PROK_ml, PROK_sem, EUK_ml, EUK_sem, FLP_L, FLP_ml),
) %>%
tab_source_note(source_note = "NAs indicate values were unavailable.
Zero values for rates indicate no grazing pressure detected.") %>%
fmt_scientific(columns = everything()) %>%
tab_options(
table.font.size = 12,
table.border.top.color = "black",
column_labels.border.bottom.color = "black",
column_labels.border.bottom.width= px(3),
table.width = pct(100))| # of incubations | Cell counts | FLPs grazer-1 | ml grazer-1 hr-1 | Specific grazing rate | Grazing rate: prokaryote cells consumed | Turnover | |||||||||||||||
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| Prokaryote cells ml-1 | SEM prokaryote cells ml-1 | Eukaryote cells ml-1 | SEM eukaryote cells ml-1 | FLP L-1 | FLP ml-1 | hour-1 | minute-1 | day-1 | ml grazer-1 hr-1 | nl grazer-1 hr-1 | Prokaryote grazer-1 hr-1 | Prokaryote grazer-1 day-1 | Prokaryote ml-1 hr-1 | MIN | MAX | Prokaryote ml-1 day-1 | Bacteria turnover % day-1 | MIN | MAX | ||
| Von Damm - Bag - Background | |||||||||||||||||||||
| Background | 3.00 | 3.79 × 104 | 8.61 × 103 | 9.18 × 101 | 2.19 × 101 | 4.86 × 106 | 4.86 × 103 | 1.78 × 10−1 | 2.96 × 10−3 | 4.26 | 3.65 × 10−5 | 3.65 × 10−11 | 1.38 | 3.32 × 101 | 1.27 × 102 | 9.68 × 101 | 1.57 × 102 | 3.05 × 103 | 8.05 | 1.04 × 101 | 6.56 |
| Plume | 3.00 | 1.65 × 104 | 2.62 × 103 | 1.58 × 102 | 6.71 × 101 | 3.42 × 107 | 3.42 × 104 | 3.16 × 10−1 | 5.27 × 10−3 | 7.59 | 9.24 × 10−6 | 9.24 × 10−12 | 1.52 × 10−1 | 3.65 | 2.40 × 101 | 1.38 × 101 | 3.42 × 101 | 5.77 × 102 | 3.50 | 4.16 | 3.02 |
| Von Damm - Bag - Vent fluid | |||||||||||||||||||||
| X-18 | 2.00 | 1.11 × 105 | 2.97 × 103 | 3.15 × 102 | 1.05 × 102 | 3.15 × 106 | 3.15 × 103 | 1.05 × 10−1 | 1.74 × 10−3 | 2.51 | 3.32 × 10−5 | 3.32 × 10−11 | 3.70 | 8.89 × 101 | 1.17 × 103 | 7.78 × 102 | 1.56 × 103 | 2.80 × 104 | 2.51 × 101 | 2.58 × 101 | 2.45 × 101 |
| Ravelin #2 | 3.00 | 7.11 × 104 | NA | 4.09 × 102 | 7.35 × 101 | 5.19 × 107 | 5.19 × 104 | 2.08 × 10−1 | 3.47 × 10−3 | 5.00 | 4.01 × 10−6 | 4.01 × 10−12 | 2.85 × 10−1 | 6.85 | 1.17 × 102 | 9.59 × 101 | 1.38 × 102 | 2.80 × 103 | 3.94 | NA | NA |
| Mustard Stand | 2.00 | 5.67 × 104 | 1.43 × 104 | 2.60 × 102 | 2.89 × 101 | 6.60 × 106 | 6.60 × 103 | 0.00 | 0.00 | 0.00 | 0.00 | 0.00 | 0.00 | 0.00 | 0.00 | 0.00 | 0.00 | 0.00 | 0.00 | 0.00 | 0.00 |
| Shrimp Hole | 2.00 | 4.20 × 104 | 3.15 × 103 | 3.86 × 102 | 7.87 | 1.13 × 108 | 1.13 × 105 | 0.00 | 0.00 | 0.00 | 0.00 | 0.00 | 0.00 | 0.00 | 0.00 | 0.00 | 0.00 | 0.00 | 0.00 | 0.00 | 0.00 |
| Von Damm - IGT - Vent fluid | |||||||||||||||||||||
| Old Man Tree | 2.00 | 7.11 × 104 | NA | 4.72 × 102 | 1.22 × 102 | 5.35 × 106 | 5.35 × 103 | 8.71 × 10−1 | 1.45 × 10−2 | 2.09 × 101 | 1.63 × 10−4 | 1.63 × 10−10 | 1.16 × 101 | 2.78 × 102 | 5.47 × 103 | 4.05 × 103 | 6.88 × 103 | 1.31 × 105 | 1.84 × 102 | NA | NA |
| Ravelin #2 | 2.00 | 7.11 × 104 | NA | 6.21 × 102 | 1.23 × 102 | 5.35 × 106 | 5.35 × 103 | 9.24 × 10−1 | 1.54 × 10−2 | 2.22 × 101 | 1.73 × 10−4 | 1.73 × 10−10 | 1.23 × 101 | 2.95 × 102 | 7.62 × 103 | 6.11 × 103 | 9.14 × 103 | 1.83 × 105 | 2.57 × 102 | NA | NA |
| Ravelin #2 | 2.00 | 7.11 × 104 | NA | 6.21 × 102 | 1.23 × 102 | 5.35 × 106 | 5.35 × 103 | 9.62 × 10−2 | 1.60 × 10−3 | 2.31 | 1.80 × 10−5 | 1.80 × 10−11 | 1.28 | 3.07 × 101 | 7.94 × 102 | 6.36 × 102 | 9.51 × 102 | 1.91 × 104 | 2.68 × 101 | NA | NA |
| Piccard - Bag - Background | |||||||||||||||||||||
| Plume | 3.00 | 5.14 × 104 | 4.62 × 103 | 7.93 × 101 | 1.68 × 101 | 2.97 × 107 | 2.97 × 104 | 3.22 × 10−1 | 5.36 × 10−3 | 7.72 | 1.09 × 10−5 | 1.09 × 10−11 | 5.58 × 10−1 | 1.34 × 101 | 4.43 × 101 | 3.49 × 101 | 5.36 × 101 | 1.06 × 103 | 2.07 | 2.27 | 1.90 |
| Piccard - Bag - Vent fluid | |||||||||||||||||||||
| Shrimpocalypse | 1.00 | 2.39 × 105 | 6.58 × 104 | 4.55 × 102 | NA | 1.70 × 107 | 1.70 × 104 | 9.41 × 10−1 | 1.57 × 10−2 | 2.26 × 101 | 5.54 × 10−5 | 5.54 × 10−11 | 1.32 × 101 | 3.17 × 102 | 6.01 × 103 | NA | NA | 1.44 × 105 | 6.04 × 101 | 8.34 × 101 | 4.74 × 101 |
| Lots 'O Shrimp | 1.00 | 5.39 × 104 | 1.37 × 104 | 2.31 × 102 | NA | 6.17 × 105 | 6.17 × 102 | 0.00 | 0.00 | 0.00 | 0.00 | 0.00 | 0.00 | 0.00 | 0.00 | NA | NA | 0.00 | 0.00 | 0.00 | 0.00 |
| Piccard - IGT - Vent fluid | |||||||||||||||||||||
| Shrimpocalypse | 1.00 | 2.39 × 105 | 6.58 × 104 | 4.55 × 102 | 7.00 × 101 | 5.35 × 106 | 5.35 × 103 | 1.01 | 1.68 × 10−2 | 2.42 × 101 | 1.88 × 10−4 | 1.88 × 10−10 | 4.49 × 101 | 1.08 × 103 | 2.04 × 104 | 1.73 × 104 | 2.36 × 104 | 4.90 × 105 | 2.05 × 102 | 2.84 × 102 | 1.61 × 102 |
| NAs indicate values were unavailable. Zero values for rates indicate no grazing pressure detected. | |||||||||||||||||||||
8 Estimates of carbon biomass and consumption
8.1 Protist carbon-measured
References for estimating biovolume Pernice, M.C., Forn, I., Gomes, A., Lara, E., Alonso-Sáez, L., Arrieta, J.M., et al. (2015) Global abundance of planktonic heterotrophic protists in the deep ocean. ISME J 9: 782–792.
# Import manual biovolume measurements
biov <- read.delim("input-data/biovol-euk-12-10-2020.txt")
# head(biov)
# Calculate volume
biov_calc <- biov %>%
mutate(SizeFrac = case_when(
h >= 20 ~ "micro",
TRUE ~ "nano")) %>%
mutate(Volume = ((pi/6) * (d^2) * d)) %>% # Calculate volume (um cubed) # Hillebrand et al. 1999
mutate(pgC_cell = (0.216 * (Volume^0.939))) %>% # Calculate Cell biomass in pg C per cell # Menden-Deuer and Lessard 2000
data.frame
# View(biov_calc)
biov_calc## EXP VENT_BSW h d SizeFrac Volume pgC_cell
## 1 IGT vent 30.077 25.764 micro 8954.44130 1110.2426245
## 2 IGT vent 89.582 10.000 micro 523.59878 77.1957618
## 3 Bag BSW 14.595 8.036 nano 271.71800 41.6956679
## 4 Bag BSW 12.480 8.982 nano 379.41786 57.0486460
## 5 Bag vent 9.218 3.120 nano 15.90239 2.9015292
## 6 IGT vent 17.255 9.986 nano 521.40274 76.8917043
## 7 Bag vent 41.153 21.000 micro 4849.04826 624.1445904
## 8 IGT vent 10.282 4.136 nano 37.04591 6.4194942
## 9 IGT vent 29.776 25.852 micro 9046.50993 1120.9583343
## 10 IGT vent 10.991 4.000 nano 33.51032 5.8424695
## 11 Bag vent 14.333 2.000 nano 4.18879 0.8290772
## 12 Bag vent 36.164 3.000 micro 14.13717 2.5980292
## 13 Bag BSW 16.206 14.924 nano 1740.42111 238.4669404
## 14 Bag BSW 7.000 7.000 nano 179.59438 28.2640658
## 15 Bag vent 10.069 5.000 nano 65.44985 10.9544849Volume is reported as um^3
# Volume by experiment type
biov_calc %>%
group_by(EXP) %>% summarise(VOL = mean(Volume), C = mean(pgC_cell))## # A tibble: 2 × 3
## EXP VOL C
## <chr> <dbl> <dbl>
## 1 Bag 836. 112.
## 2 IGT 3186. 400.# Volume by euk size
biov_calc %>%
group_by(SizeFrac) %>% summarise(VOL = mean(Volume), C = mean(pgC_cell))## # A tibble: 2 × 3
## SizeFrac VOL C
## <chr> <dbl> <dbl>
## 1 micro 4678. 587.
## 2 nano 325. 46.9# Volume by site
biov_calc %>%
group_by(VENT_BSW) %>% summarise(VOL = mean(Volume), C = mean(pgC_cell))## # A tibble: 2 × 3
## VENT_BSW VOL C
## <chr> <dbl> <dbl>
## 1 BSW 643. 91.4
## 2 vent 2188. 276.# head(biov_calc)
euk_vol <- mean(biov_calc$Volume);euk_vol # in um^3## [1] 1775.759euk_carbon <- mean(biov_calc$pgC_cell); euk_carbon # in pg C per cell## [1] 226.9636euk_carbon_min <- min(biov_calc$pgC_cell); euk_carbon_min## [1] 0.8290772euk_carbon_max <- max(biov_calc$pgC_cell); euk_carbon_max## [1] 1120.958# euk_carbonAvg euk biomass pg C per individual cell == {r}euk_carbon
8.2 Protist carbon-literature
Compare with Menden-Deuer and Lessard 2000, Table 2 - using only the heterotrophic species measured. Based on Table 2, the min volume was 4745 and the maximum was 1.2 x10^7 µm^3. Carbon content was measured at pg per cell, this was 469.48-35,339 pg per cell.
Import the heterotroph species volume and carbon content to compare to my measured values.
# Hu-measured
range(biov_calc$Volume)## [1] 4.18879 9046.50993range(biov_calc$pgC_cell)## [1] 0.8290772 1120.9583343c_prev <- read.delim("input-data/md-lessard-2000.txt") # Table 2, heterotrophs only
# combine and plot
carbon_compare <- c_prev %>%
add_column(source = "Menden-Deuer Lessard") %>%
select(source, Volume = vol, pgC_cell) %>%
rbind(biov_calc %>% add_column(source = "MCR") %>% select(source, Volume, pgC_cell)) %>%
ggplot(aes(x = Volume, y = pgC_cell, fill = source)) +
geom_point(aes(fill = source), shape = 23, color = "black", size = 3) +
scale_y_log10() + scale_x_log10() +
labs(title = "Compare literature to measured cell volume & C content",
x = bquote("Volume" ~µm^-3),
y = bquote("pg C" ~cell^-1)) +
theme_bw() + theme(legend.title = element_blank(),
axis.title = element_text(size = 14),
axis.text = element_text(size = 14),
legend.text = element_text(size = 14))
carbon_compare
Upon comparison, the measured carbon content was much lower from the grazing experiments. This makes sense, as I am looking at preserved specimen and a smaller total number of cells. AND the deep-sea protist cell sizes may be smaller overall.
Find lowest estimates or protist carbon, benthic estimates, and others? How does it compare to my measurements?
euk_carbon_lit_mean <- mean(c_prev$pgC_cell)
euk_carbon_lit_min <- min(c_prev$pgC_cell)
euk_carbon_lit_max <- max(c_prev$pgC_cell)8.3 Prokaryote carbon-literature
Below adding in biomass estimates from prokaryotes and protists.
bac_carbon_ug <- (86)*(1.00E-9) # From Derived from Morono et al. 2011
# bac_carbon_ug
bac_carbon_ug_2 <- (173)*(1.00E-9) # Derived from McNichol et al. 2018; LOFERER-KRO ̈ ßBACHER, J. KLIMA & R. PSENNER 1998
# table_wcalcs8.4 Generate carbon biomass table
Incorporate calculations that include biomass of population and ug C consumed. For rate measurements, only incorporate the Morono et al. 2011 biomass for prokaryotes. This way it is on the lower end and is comparable to Gorda Ridge work.
NEED TO MODIFY SO IT INCLUDES BIOMASS OF EUK AND PROK AT EACH SITE
bsw <- c("Plume", "Background")
table_wcalcs_biomass <- table_wcalcs %>%
add_column(euk_C_ug_Hu = (euk_carbon / (1.00E+06))) %>% # Convert to ug from pg
add_column(euk_C_ug_lit = (euk_carbon_lit_mean / (1.00E+06))) %>% # literature
add_column(bac_C_ug = bac_carbon_ug) %>%
add_column(bac_C_ug_2 = bac_carbon_ug_2) %>%
# Grazing rate in ug C per bac per day
mutate(RATE_ugCbac_pergrazer_perday = (RATE_hr * 24 * bac_C_ug), # Grazing rate as ug C per grazer per day
# % of cell carbon per day
SPEC_INGESTION_RATE = (RATE_ugCbac_pergrazer_perday / euk_C_ug_Hu),
SPEC_INGESTION_RATE_lit = (RATE_ugCbac_pergrazer_perday / euk_C_ug_lit),
Prok_biomass = PROK_ml * bac_carbon_ug,
Euk_biomass_Hu = EUK_ml * euk_C_ug_Hu,
Euk_biomass_lit = EUK_ml * euk_C_ug_lit,
Prok_biomass_L = PROK_ml * bac_carbon_ug * 1000,
Euk_biomass_Hu_L = EUK_ml * euk_C_ug_Hu * 1000,
Euk_biomass_lit_L = EUK_ml * euk_C_ug_lit * 1000,
# Repeat with SEM values
Prok_biomass_sem = PROK_sem * bac_carbon_ug,
Euk_biomass_Hu_sem = EUK_sem * euk_C_ug_Hu,
Euk_biomass_lit_sem = EUK_sem * euk_C_ug_lit,
Prok_biomass_sem_L = PROK_sem * (bac_carbon_ug* 1000),
Euk_biomass_Hu_sem_L = EUK_sem * (euk_C_ug_Hu * 1000),
Euk_biomass_lit_sem_L = EUK_sem * (euk_C_ug_lit * 1000)) %>%
type.convert(as.is = TRUE) %>%
mutate(detected = case_when(
RATE_min < 0 ~ "Not detected",
TRUE ~ "Detected")) %>%
mutate(type = case_when(
NAME %in% bsw ~ NAME,
TRUE ~ paste("Vent", EXP, sep="-")
)) %>%
mutate(GRAZE_RATE = case_when(
RATE_min < 0 ~ 0,
TRUE ~ RATE_min
)) %>%
mutate(type_site = case_when(
NAME %in% bsw ~ NAME,
TRUE ~ "Vent"
)) %>%
data.frame
# View(table_wcalcs_biomass)- Grazing rate column == FLP per minute consumed
- Grazing effect hr == cells per ml per hr
Also make a “bounded” table that demonstrates the ug C consumed in the context of McNichol et al.
# G = number of cells grazed during experiment duration
table_wcalcs_biomass_bounded <- table_wcalcs_biomass %>%
add_column(fgC_cell = 86) %>% # Add in Morono et al. 2011 value
mutate(
# cells_consumed_perday = (G / 1), # Rate of cells consumed * in situ prok, per day
fgC_ml_perday = (GRAZING_EFFECT_day * fgC_cell), # Convert cell amount to fg C
ugC_L_perday = (fgC_ml_perday * (1e-09) * 1000), # Convert to ug C per L
lower_mcnichol = 100*(ugC_L_perday / 17.3),
upper_mcnichol = 100*(ugC_L_perday / 321.4)
) %>%
data.framehead(table_wcalcs_biomass_bounded)## SITE NAME EXP SAMPLE RATE_min REPs FLP_ml PROK_ml
## 1 Von Damm Background Bag VD-BSW 0.002958889 3 4861.824 37889.62
## 2 Von Damm Plume Bag VD-Plume 0.005274231 3 34242.354 16478.31
## 3 Von Damm X-18 Bag VD-X18 0.001744429 2 3148.722 111429.78
## 4 Von Damm Old Man Tree IGT VD-OMT 0.014510943 2 5352.828 71147.40
## 5 Von Damm Ravelin #2 Bag VD-Rav2 0.003470217 3 51890.942 71147.40
## 6 Von Damm Ravelin #2 IGT VD-Rav2 0.015395240 2 5352.828 71147.40
## PROK_sem EUK_ml EUK_sem RATE_hr RATE_day FLP_L CLEARANCE_RATE_ml
## 1 8608.427 91.83773 21.86613 0.1775333 4.260800 4861824 3.651579e-05
## 2 2623.935 157.77468 67.09859 0.3164538 7.594892 34242354 9.241591e-06
## 3 2973.793 314.87222 104.95741 0.1046657 2.511977 3148722 3.324070e-05
## 4 NA 472.30833 122.45031 0.8706566 20.895758 5352828 1.626536e-04
## 5 NA 409.33389 73.47019 0.2082130 4.997113 51890942 4.012512e-06
## 6 NA 620.99799 123.17702 0.9237144 22.169145 5352828 1.725657e-04
## CLEARANCE_RATE_nL SPEC_GRAZE_RATE_hr GRAZE_RATE_DAY GRAZING_EFFECT_hr
## 1 3.651579e-11 1.3835695 33.205668 127.06388
## 2 9.241591e-12 0.1522858 3.654860 24.02685
## 3 3.324070e-11 3.7040034 88.896081 1166.28778
## 4 1.626536e-10 11.5723785 277.737084 5465.73080
## 5 4.012512e-12 0.2854798 6.851515 116.85656
## 6 1.725657e-10 12.2775993 294.662382 7624.36451
## GRAZING_EFFECT_hr_min GRAZING_EFFECT_hr_max GRAZING_EFFECT_day
## 1 96.81058 157.31719 3049.5332
## 2 13.80868 34.24501 576.6444
## 3 777.52519 1555.05037 27990.9067
## 4 4048.68948 6882.77212 131177.5392
## 5 95.88230 137.83081 2804.5574
## 6 6112.04639 9136.68264 182984.7483
## BAC_TURNOVER_PERC BAC_TURNOVER_PERC_min BAC_TURNOVER_PERC_max euk_C_ug_Hu
## 1 8.048465 10.414647 6.558411 0.0002269636
## 2 3.499414 4.162182 3.018725 0.0002269636
## 3 25.119772 25.808540 24.466811 0.0002269636
## 4 184.374334 NA NA 0.0002269636
## 5 3.941897 NA NA 0.0002269636
## 6 257.191065 NA NA 0.0002269636
## euk_C_ug_lit bac_C_ug bac_C_ug_2 RATE_ugCbac_pergrazer_perday
## 1 0.01147298 8.6e-08 1.73e-07 3.664288e-07
## 2 0.01147298 8.6e-08 1.73e-07 6.531607e-07
## 3 0.01147298 8.6e-08 1.73e-07 2.160300e-07
## 4 0.01147298 8.6e-08 1.73e-07 1.797035e-06
## 5 0.01147298 8.6e-08 1.73e-07 4.297517e-07
## 6 0.01147298 8.6e-08 1.73e-07 1.906547e-06
## SPEC_INGESTION_RATE SPEC_INGESTION_RATE_lit Prok_biomass Euk_biomass_Hu
## 1 0.001614483 3.193840e-05 0.003258508 0.02084382
## 2 0.002877822 5.693033e-05 0.001417135 0.03580910
## 3 0.000951827 1.882946e-05 0.009582961 0.07146452
## 4 0.007917726 1.566319e-04 0.006118676 0.10719678
## 5 0.001893483 3.745771e-05 0.006118676 0.09290388
## 6 0.008400232 1.661770e-04 0.006118676 0.14094392
## Euk_biomass_lit Prok_biomass_L Euk_biomass_Hu_L Euk_biomass_lit_L
## 1 1.053653 3.258508 20.84382 1053.653
## 2 1.810146 1.417135 35.80910 1810.146
## 3 3.612524 9.582961 71.46452 3612.524
## 4 5.418786 6.118676 107.19678 5418.786
## 5 4.696281 6.118676 92.90388 4696.281
## 6 7.124700 6.118676 140.94392 7124.700
## Prok_biomass_sem Euk_biomass_Hu_sem Euk_biomass_lit_sem Prok_biomass_sem_L
## 1 0.0007403247 0.004962814 0.2508697 0.7403247
## 2 0.0002256584 0.015228935 0.7698210 0.2256584
## 3 0.0002557462 0.023821507 1.2041746 0.2557462
## 4 NA 0.027791758 1.4048704 NA
## 5 NA 0.016675055 0.8429222 NA
## 6 NA 0.027956696 1.4132080 NA
## Euk_biomass_Hu_sem_L Euk_biomass_lit_sem_L detected type GRAZE_RATE
## 1 4.962814 250.8697 Detected Background 0.002958889
## 2 15.228935 769.8210 Detected Plume 0.005274231
## 3 23.821507 1204.1746 Detected Vent-Bag 0.001744429
## 4 27.791758 1404.8704 Detected Vent-IGT 0.014510943
## 5 16.675055 842.9222 Detected Vent-Bag 0.003470217
## 6 27.956696 1413.2080 Detected Vent-IGT 0.015395240
## type_site fgC_cell fgC_ml_perday ugC_L_perday lower_mcnichol upper_mcnichol
## 1 Background 86 262259.86 0.26225986 1.5159529 0.08159921
## 2 Plume 86 49591.42 0.04959142 0.2866556 0.01542981
## 3 Vent 86 2407217.98 2.40721798 13.9145548 0.74897884
## 4 Vent 86 11281268.37 11.28126837 65.2096438 3.51003994
## 5 Vent 86 241191.93 0.24119193 1.3941730 0.07504416
## 6 Vent 86 15736688.36 15.73668836 90.9635165 4.89629383# View(table_wcalcs_biomass_bounded)
# write_delim(table_wcalcs_biomass_bounded, file = "table-wcalc.txt", delim = "\t")
# ?write_delim()8.5 Compilation of plots
Format and factor values to plot. * Grazing rate column == FLP per minute consumed * Grazing effect hr == cells per ml per hr * Specific ingestion rate == % of cell carbon per day, estimated with literature and measured carbon values
biomass_rate_plot <- table_wcalcs_biomass_bounded %>%
select(SITE, NAME, EXP, SAMPLE, type,
PROK_ml, EUK_ml, PROK_sem, EUK_sem,
Prok_biomass_L, Euk_biomass_Hu_L, Euk_biomass_lit_L,
Prok_biomass_sem_L, Euk_biomass_Hu_sem_L, Euk_biomass_lit_sem_L,
GRAZING_EFFECT_hr, GRAZING_EFFECT_hr_min, GRAZING_EFFECT_hr_max, BAC_TURNOVER_PERC, BAC_TURNOVER_PERC_min, BAC_TURNOVER_PERC_max,
GRAZE_RATE, RATE_ugCbac_pergrazer_perday,
SPEC_INGESTION_RATE, SPEC_INGESTION_RATE_lit) %>%
pivot_longer(cols = c(PROK_ml, EUK_ml,
Prok_biomass_L, Euk_biomass_Hu_L, Euk_biomass_lit_L,
GRAZING_EFFECT_hr, BAC_TURNOVER_PERC,
GRAZE_RATE, RATE_ugCbac_pergrazer_perday,
SPEC_INGESTION_RATE, SPEC_INGESTION_RATE_lit),
names_to = "Variable", values_to = "Value") %>%
pivot_longer(cols = c(PROK_sem, EUK_sem, Prok_biomass_sem_L, Euk_biomass_Hu_sem_L, Euk_biomass_lit_sem_L,
BAC_TURNOVER_PERC_max, BAC_TURNOVER_PERC_min, GRAZING_EFFECT_hr_max, GRAZING_EFFECT_hr_min, GRAZING_EFFECT_hr_max),
names_to = "SEM_variable", values_to = "SEM") %>%
data.frame
biomass_rate_plot$NAME_ORDER <- factor(biomass_rate_plot$NAME, levels = c("Background","Plume", "Quakeplume", "Shrimpocalypse", "Lots 'O Shrimp", "X-18", "Old Man Tree", "Ravelin #2", "Mustard Stand", "Shrimp Hole", "Hot Chimlet #1", "Shrimp Gulley", "South of Hot Chimlet", "South of LungSnack", "Arrow Loop", "Bartizan", "Ravelin #1"))
biomass_rate_plot$VARIABLE_ORDER <- factor(biomass_rate_plot$Variable,
levels = c("PROK_ml", "EUK_ml",
"Prok_biomass_L", "Euk_biomass_Hu_L", "Euk_biomass_lit_L",
"GRAZING_EFFECT_hr", "BAC_TURNOVER_PERC",
"GRAZE_RATE", "RATE_ugCbac_pergrazer_perday",
"SPEC_INGESTION_RATE", "SPEC_INGESTION_RATE_lit"),
labels = c("Prokaryote~cells~mL^{-1}", "Eukaryote~cells~mL^{-1}",
"Prokaryote~µg~C~L^{-1}", "Measured~eukaryote~µg~C~L^{-1}", "Literature-based~eukaryote~µg~C~L^{-1}",
"Cells~mL^{-1}~hr^{-1}", "Prokaryote~turnover~'%'~d^{-1}",
"FLP~consumed~min^{-1}", "µg~C~grazer^{-1}~day^{-1}",
"Cell~carbon~%~day^{-1}", "Cell~carbon~%~day^{-1}~(lit)"))
# View(biomass_rate_plot)Function to generate consistent plots for Mid-Cayman Rise cell concentration and biomass data.
conc_rate_plot_mcr <- function(df, var, sem){
df %>%
filter(Variable == var) %>%
filter(SEM_variable == sem) %>%
ggplot(aes(y = Value, x = NAME_ORDER, shape = EXP, fill = SITE)) +
geom_errorbar(aes(ymax = (Value + SEM), ymin = (Value - SEM)),
width = 0.2, position = position_dodge(width = 0.4)) +
geom_point(stat = "identity", aes(shape = EXP, fill = SITE),
color = "black", size = 3, position = position_dodge(width = 0.4)) +
scale_shape_manual(values = c(21, 23)) +
scale_fill_manual(values = c("#de2d26", "#1c9099")) +
facet_wrap(VARIABLE_ORDER ~ ., scales = "free",
strip.position = c("left"), labeller = label_parsed) +
scale_y_log10() +
# scale_y_log10(labels = function(x) format(x, scientific = TRUE)) +
theme_minimal() +
theme(panel.grid.major = element_line(), panel.grid.minor = element_blank(),
panel.background = element_blank(),
axis.line = element_line(colour = "black"),
axis.text.x = element_text(color="black", size = 11,
angle = 45, hjust = 1, vjust = 1),
axis.text.y = element_text(color="black", size = 11),
axis.title =element_text(color="black", size = 14),
axis.ticks = element_line(),
legend.title = element_blank(),
strip.placement = "outside",
strip.text.y = element_text(color="black", size = 11),
strip.text.x = element_blank())+
guides(fill = guide_legend(override.aes = list(shape = c(21))),
shape = guide_legend(override.aes = list(fill = "black"))) +
labs(x = "", y = "")
}# svg("cells-per-ml-yaxies.svg", h = 8, w = 13)
plot_grid(conc_rate_plot_mcr(biomass_rate_plot, "PROK_ml", "PROK_sem"),
conc_rate_plot_mcr(biomass_rate_plot, "EUK_ml", "EUK_sem"),
conc_rate_plot_mcr(biomass_rate_plot, "Prok_biomass_L", "Prok_biomass_sem_L"),
conc_rate_plot_mcr(biomass_rate_plot, "Euk_biomass_Hu_L", "Euk_biomass_Hu_sem_L"))
# dev.off()
# svg("MCR-only-cellcounts-grazingrate.svg", h=12, w=7)
plot_grid(conc_rate_plot_mcr(biomass_rate_plot, "PROK_ml", "PROK_sem"),
conc_rate_plot_mcr(biomass_rate_plot, "EUK_ml", "EUK_sem"),
conc_rate_plot_mcr(biomass_rate_plot, "GRAZING_EFFECT_hr", "Prok_biomass_sem_L"), ncol = 1)## Warning: Transformation introduced infinite values in continuous y-axis## Warning in self$trans$transform(x): NaNs produced## Warning: Transformation introduced infinite values in continuous y-axis
## Warning: Transformation introduced infinite values in continuous y-axis
# dev.off()subset <- c("GRAZING_EFFECT_hr", "RATE_ugCbac_pergrazer_perday")
biomass_rate_plot %>%
filter(Variable %in% subset) %>%
group_by(EXP, type, Variable) %>%
summarise(MEAN = mean(Value))## `summarise()` has grouped output by 'EXP', 'type'. You can override using the `.groups` argument.## # A tibble: 8 × 4
## # Groups: EXP, type [4]
## EXP type Variable MEAN
## <chr> <chr> <chr> <dbl>
## 1 Bag Background GRAZING_EFFECT_hr 1.27e+2
## 2 Bag Background RATE_ugCbac_pergrazer_perday 3.66e-7
## 3 Bag Plume GRAZING_EFFECT_hr 3.41e+1
## 4 Bag Plume RATE_ugCbac_pergrazer_perday 6.59e-7
## 5 Bag Vent-Bag GRAZING_EFFECT_hr 1.21e+3
## 6 Bag Vent-Bag RATE_ugCbac_pergrazer_perday 4.31e-7
## 7 IGT Vent-IGT GRAZING_EFFECT_hr 8.58e+3
## 8 IGT Vent-IGT RATE_ugCbac_pergrazer_perday 1.50e-6# 8.577791e+03 - 1.214981e+03 # scale_y_continuous(labels = function(x) format(x, scientific = TRUE)) +
grazing_rate <- biomass_rate_plot %>%
type.convert(as.is = TRUE) %>%
filter(Variable == "GRAZING_EFFECT_hr") %>%
filter(SEM_variable == "GRAZING_EFFECT_hr_min" | SEM_variable == "GRAZING_EFFECT_hr_max") %>%
pivot_wider(names_from = SEM_variable, values_from = SEM) %>%
ggplot(aes(y = Value, x = NAME_ORDER, shape = EXP, fill = SITE)) +
geom_errorbar(aes(ymax = (GRAZING_EFFECT_hr_max), ymin = (GRAZING_EFFECT_hr_min)),
width = 0.2, position = position_dodge(width = 0.4)) +
geom_point(stat = "identity", aes(shape = EXP, fill = SITE),
color = "black", size = 3, position = position_dodge(width = 0.4)) +
scale_shape_manual(values = c(21, 24)) +
scale_fill_manual(values = c("#de2d26", "#1c9099")) +
facet_wrap(VARIABLE_ORDER ~ ., scales = "free",
strip.position = c("left"), labeller = label_parsed) +
theme_minimal() +
theme(panel.grid.major = element_line(), panel.grid.minor = element_blank(),
panel.background = element_blank(),
axis.line = element_line(colour = "black"),
axis.text.x = element_text(color="black", size = 11,
angle = 45, hjust = 1, vjust = 1),
axis.text.y = element_text(color="black", size = 11),
axis.title =element_text(color="black", size = 11),
axis.ticks = element_line(),
legend.title = element_blank(),
strip.placement = "outside",
strip.text.y = element_text(color="black", size = 11),
strip.text.x = element_blank())+
guides(fill = guide_legend(override.aes = list(shape = c(21))),
shape = guide_legend(override.aes = list(fill = "black"))) +
labs(x = "", y = "")
# class(grazing_rate)
# head(biomass_rate_plot)
# svg("plot-grazing-withscale.svg", h = 7, w = 6)
# grazing_rate
plot_grid(grazing_rate + scale_y_continuous(limits = c(0,2000)),
grazing_rate,
ncol = 1, rel_heights = c(0.6,1))## Warning: Removed 4 rows containing missing values (geom_point).
# dev.off()# scale_y_continuous(labels = function(x) format(x, scientific = TRUE)) +
prok_turnover <- biomass_rate_plot %>%
filter(Variable == "BAC_TURNOVER_PERC") %>%
filter(SEM_variable == "BAC_TURNOVER_PERC_min" | SEM_variable == "BAC_TURNOVER_PERC_max") %>%
pivot_wider(names_from = SEM_variable, values_from = SEM) %>%
ggplot(aes(y = Value, x = NAME_ORDER, fill = SITE)) +
geom_bar(stat = "identity", aes(fill = SITE),
color = "black", width = 2, position = position_dodge2(preserve = "single")) +
geom_errorbar(aes(ymax = (BAC_TURNOVER_PERC_max), ymin = (BAC_TURNOVER_PERC_min)),
width = 0.3, position = position_dodge2(preserve = "single")) +
scale_fill_manual(values = c("#de2d26", "#1c9099")) +
facet_wrap(VARIABLE_ORDER ~ ., scales = "free",
strip.position = c("left"), labeller = label_parsed) +
theme_minimal() +
theme(panel.grid.major = element_line(), panel.grid.minor = element_blank(),
panel.background = element_blank(),
axis.line = element_line(colour = "black"),
axis.text.x = element_text(color="black", size = 11,
angle = 45, hjust = 1, vjust = 1),
axis.text.y = element_text(color="black", size = 11),
axis.title =element_text(color="black", size = 11),
axis.ticks = element_line(),
legend.title = element_blank(),
strip.placement = "outside",
strip.text.y = element_text(color="black", size = 11),
strip.text.x = element_blank())+
guides(fill = guide_legend(override.aes = list(shape = c(21))),
shape = guide_legend(override.aes = list(fill = "black"))) +
labs(x = "", y = "")# View(biomass_rate_plot)
# tmp <- (biomass_rate_plot %>% filter(Variable == "RATE_ugCbac_pergrazer_perday" & Value > 0))
# range(tmp$Value)
# View(tmp)
# unique(biomass_rate_plot$Variable)# svg("prok-turnover.svg", h = 8, w = 6)
plot_grid(prok_turnover,
prok_turnover + scale_y_continuous(limits = c(0,100)),
ncol = 1)## Warning: Removed 3 rows containing missing values (geom_bar).
# dev.off()8.6 Compare MCR data with Gorda Ridge - think about environmental parameters
Import Gorda Ridge data
gr <- read.delim("../../GordaRidgeCruise_2019/protist-gordaridge-2021/Grazing-at-GordaRidge-SKH-2021/data-input/Grazing-calc-wCarbon-results.txt")
# env_gr <- read.delim("../../GordaRidgeCruise_2019/protist-gordaridge-2021/Grazing-at-GordaRidge-SKH-2021/data-input/GR-environ-SAMPLE.txt")
temps <- read.delim("/Users/sarahhu/Desktop/Projects/GordaRidgeCruise_2019/protist-gordaridge-2019-analysis/temperature-allvents.txt")
mcr_graze <- read.delim("table-wcalc.txt")
# head(mcr_graze)all_vents <- mcr_graze %>%
type.convert(as.is = TRUE) %>%
select(SITE, NAME, SAMPLE, EXP, PROK_ml, EUK_ml, GRAZING_EFFECT_hr, GRAZING_EFFECT_hr_min, GRAZING_EFFECT_hr_max, BAC_TURNOVER_PERC, ugC_L_perday) %>%
rbind(gr %>%
add_column(SITE = "Gorda Ridge") %>%
add_column(EUK_ml = NA) %>%
separate(SAMPLE, c("SAMPLE", "NAME"), sep = "-") %>%
select(SITE, NAME, SAMPLE, EXP = Bottle, PROK_ml = prok_avg, EUK_ml, GRAZING_EFFECT_hr = GrazingRate_hr, GRAZING_EFFECT_hr_min = GrazingRate_hr_min, GRAZING_EFFECT_hr_max = GrazingRate_hr_max, BAC_TURNOVER_PERC = Prok_turnover, ugC_L_perday = ugC_L_perday_morono)) %>%
left_join(temps) %>%
mutate(SAMPLE_TYPE = case_when(
grepl("BSW", NAME) ~ "Background",
grepl("Near vent BW", NAME) ~ "Background",
grepl("Background", NAME) ~ "Background",
grepl("Plume", NAME) ~ "Background",
TRUE ~ "Vent"
))## Joining, by = c("NAME", "SAMPLE")# View(all_vents)
# all_vents
# write_delim(all_vents, file = "grazing-cellcounts-GR_MCR.txt", delim = "\t")Plot grazing rates
allrates <- all_vents %>%
select(SITE, NAME, SAMPLE, EXP, SAMPLE_TYPE, starts_with("GRAZING_EFFECT_")) %>%
distinct() %>%
ggplot(aes(y = GRAZING_EFFECT_hr, x = NAME, fill = SITE, shape = SAMPLE_TYPE)) +
geom_errorbar(aes(ymax = (GRAZING_EFFECT_hr_max), ymin = (GRAZING_EFFECT_hr_min)),
width = 0.2, position = position_dodge(width = 0.4)) +
geom_point(stat = "identity", aes(fill = SITE, shape = SAMPLE_TYPE),
color = "black", size = 3, position = position_dodge(width = 0.4)) +
scale_shape_manual(values = c(21, 24)) +
scale_fill_manual(values = c("#de2d26", "#1c9099", "#addd8e")) +
facet_grid(. ~ SAMPLE_TYPE, scales = "free", space = "free") +
theme_minimal() +
theme(panel.grid.major = element_line(), panel.grid.minor = element_blank(),
panel.background = element_blank(),
axis.line = element_line(colour = "black"),
axis.text.x = element_text(color="black", size = 11,
angle = 45, hjust = 1, vjust = 1),
axis.text.y = element_text(color="black", size = 11),
axis.title =element_text(color="black", size = 11),
axis.ticks = element_line(),
legend.title = element_blank(),
strip.placement = "outside",
strip.text.y = element_text(color="black", size = 11),
strip.text.x = element_blank())+
guides(fill = guide_legend(override.aes = list(shape = c(21))),
shape = guide_legend(override.aes = list(fill = "black"))) +
labs(x = "", y = bquote("cells"~ml^-1~hr^-1))# svg("compare-all-rates.svg", h = 10, w = 7)
plot_grid(allrates,
allrates + scale_y_continuous(limits = c(0,1000)),
allrates + scale_y_log10(), ncol = 1)## Warning: Removed 7 rows containing missing values (geom_point).## Warning: Transformation introduced infinite values in continuous y-axis
## Warning: Transformation introduced infinite values in continuous y-axis
## Warning: Transformation introduced infinite values in continuous y-axis
## Warning: Transformation introduced infinite values in continuous y-axis
# dev.off()
# svg("compare-all-rates-color.svg", h = 4, w = 9)
# allrates + scale_y_log10()
# dev.off()Repeat grazing rate plot, but removed undetectable
# unique(all_vents$GRAZING_EFFECT_hr)
allrates_nonzero <- all_vents %>%
filter(GRAZING_EFFECT_hr > 0) %>%
select(SITE, NAME, SAMPLE, EXP, SAMPLE_TYPE, starts_with("GRAZING_EFFECT_")) %>%
distinct() %>%
ggplot(aes(y = GRAZING_EFFECT_hr, x = NAME, fill = SITE, shape = SAMPLE_TYPE)) +
geom_errorbar(aes(ymax = (GRAZING_EFFECT_hr_max), ymin = (GRAZING_EFFECT_hr_min)),
width = 0.2, position = position_dodge(width = 0.4)) +
geom_point(stat = "identity", aes(fill = SITE, shape = SAMPLE_TYPE),
color = "black", size = 3, position = position_dodge(width = 0.4)) +
scale_shape_manual(values = c(21, 24)) +
scale_fill_manual(values = c("#de2d26", "#1c9099", "#addd8e")) +
facet_grid(. ~ SAMPLE_TYPE, scales = "free", space = "free") +
theme_minimal() +
theme(panel.grid.major = element_line(), panel.grid.minor = element_blank(),
panel.background = element_blank(),
axis.line = element_line(colour = "black"),
axis.text.x = element_text(color="black", size = 11,
angle = 45, hjust = 1, vjust = 1),
axis.text.y = element_text(color="black", size = 11),
axis.title =element_text(color="black", size = 11),
axis.ticks = element_line(),
legend.title = element_blank(),
strip.placement = "outside",
strip.text.y = element_text(color="black", size = 11),
strip.text.x = element_blank())+
guides(fill = guide_legend(override.aes = list(shape = c(21))),
shape = guide_legend(override.aes = list(fill = "black"))) +
labs(x = "", y = bquote("cells"~ml^-1~hr^-1))# svg("compare-all-rates-color-nonZero.svg", h = 4, w = 7)
allrates_nonzero + scale_y_log10()
# dev.off()# svg("plot-gcellconc-grazing-log.svg", h = 10, w = 6)
# grazing_rate + scale_y_log10(limits = c(1e-1,1e5), breaks=10^(0:7))
plot_grid(conc_rate_plot_mcr(biomass_rate_plot, "PROK_ml", "PROK_sem"),
conc_rate_plot_mcr(biomass_rate_plot, "EUK_ml", "EUK_sem"),
allrates_nonzero + scale_y_log10(), ncol = 1)
# dev.off()8.7 Run regression analysis with MCR and GR data
metadata <- read.delim("input-data/flp-exp-metadata-compiled.txt")
# head(metadata)
# temps
# all_ventsLinear regression with both Gorda Ridge and MCR data
library(broom)
# ?pivot_longer
# e IGT
# head(all_vents)
regression_input <- all_vents %>%
filter(!(SAMPLE == "Piccard-Shrimpocalypse" & EXP == "IGT")) %>%
filter(!is.na(Highest.Temp)) %>%
select(SITE, NAME, SAMPLE, SAMPLE_TYPE, EXP, PROK_ml, EUK_ml, GRAZING_EFFECT_hr, BAC_TURNOVER_PERC, ugC_L_perday, TEMP = Highest.Temp) %>%
mutate(PROK_EUK_RATIO = (PROK_ml/EUK_ml)) %>%
pivot_longer(cols = c(GRAZING_EFFECT_hr, BAC_TURNOVER_PERC, ugC_L_perday), names_to = "RATE", values_to = "RATE_VALUE") %>%
pivot_longer(cols = c(PROK_ml, EUK_ml, PROK_EUK_RATIO, TEMP), names_to = "PARAMS", values_to = "PARAMS_VALUE") %>%
data.frame
regression_tmp <- regression_input %>%
# Set up the linear regression
group_by(RATE, PARAMS) %>%
nest(-RATE, -PARAMS) %>%
mutate(lm_fit = map(data, ~lm(RATE_VALUE ~ PARAMS_VALUE, data = .)),
tidied = map(lm_fit, tidy)) %>%
unnest(tidied) %>%
select(RATE, PARAMS,
term, estimate) %>%
pivot_wider(names_from = term, values_from = estimate) %>%
select(everything(), SLOPE = PARAMS_VALUE) %>%
data.frame## Warning: All elements of `...` must be named.
## Did you want `data = c(SITE, NAME, SAMPLE, SAMPLE_TYPE, EXP, RATE_VALUE, PARAMS_VALUE)`?regression_results <- regression_input %>%
group_by(RATE, PARAMS) %>%
nest(-RATE, -PARAMS) %>%
mutate(lm_fit = map(data, ~lm(RATE_VALUE ~ PARAMS_VALUE, data = .)),
glanced = map(lm_fit, glance)) %>%
unnest(glanced) %>%
select(RATE, PARAMS, r.squared, adj.r.squared) %>%
right_join(regression_tmp) %>%
right_join(regression_input) %>%
data.frame## Warning: All elements of `...` must be named.
## Did you want `data = c(SITE, NAME, SAMPLE, SAMPLE_TYPE, EXP, RATE_VALUE, PARAMS_VALUE)`?## Joining, by = c("RATE", "PARAMS")
## Joining, by = c("RATE", "PARAMS")# head(regression_results)Plot regression results
regression_results %>%
# filter(RATE == "GRAZING_EFFECT_hr") %>%
# filter(PARAMS == "TEMP") %>%
ggplot(aes(x = PARAMS_VALUE, y = RATE_VALUE, shape = SAMPLE_TYPE, fill = SITE)) +
geom_abline(aes(slope = SLOPE, intercept = `X.Intercept.`), color = "black", linetype = "dashed", size = 0.5) +
geom_point(color = "black", aes(shape = SAMPLE_TYPE, fill = SITE)) +
scale_shape_manual(values = c(21, 24)) +
scale_fill_manual(values = c("#476AA7","#7299CE", "#A2937A")) +
facet_wrap(PARAMS ~ RATE + round(r.squared, 3), scales = "free", ncol = 5,
strip.position = "bottom", labeller = label_parsed) +
theme_bw() +
theme(
strip.background = element_blank(),
strip.placement = "outside",
strip.text = element_text(color = "black", size = 10),
axis.title = element_text(color = "black", size = 10),
legend.title = element_blank()) # +## Warning: Removed 24 rows containing missing values (geom_point).
# labs(y = bquote("Cells "~mL^-1 ~hr^-1), x = "Temperature (C)")NEXT STEP = get more biogeochem data?
9 Parameter comparions to plot
head(regression_results)## RATE PARAMS r.squared adj.r.squared X.Intercept. SLOPE
## 1 GRAZING_EFFECT_hr PROK_ml 0.3260524 0.2864084 -488.8238 0.02782789
## 2 GRAZING_EFFECT_hr PROK_ml 0.3260524 0.2864084 -488.8238 0.02782789
## 3 GRAZING_EFFECT_hr PROK_ml 0.3260524 0.2864084 -488.8238 0.02782789
## 4 GRAZING_EFFECT_hr PROK_ml 0.3260524 0.2864084 -488.8238 0.02782789
## 5 GRAZING_EFFECT_hr PROK_ml 0.3260524 0.2864084 -488.8238 0.02782789
## 6 GRAZING_EFFECT_hr PROK_ml 0.3260524 0.2864084 -488.8238 0.02782789
## SITE NAME SAMPLE SAMPLE_TYPE EXP RATE_VALUE PARAMS_VALUE
## 1 Von Damm Background VD-BSW Background Bag 127.06388 37889.62
## 2 Von Damm Background VD-BSW Background Bag 127.06388 37889.62
## 3 Von Damm Background VD-BSW Background Bag 127.06388 37889.62
## 4 Von Damm Plume VD-Plume Background Bag 24.02685 16478.31
## 5 Von Damm X-18 VD-X18 Vent Bag 1166.28778 111429.78
## 6 Von Damm Old Man Tree VD-OMT Vent IGT 5465.73080 71147.40unique(regression_results$PARAMS)## [1] "PROK_ml" "EUK_ml" "PROK_EUK_RATIO" "TEMP"Plot temperature by euk:prok ratio? individual? other things?
10 Prediction of MCR data with parameters
From Vaque et al. 1994 paper equation (1):
log(GT) = -3.21 + 0.99 log (HNF) + 0.028 (T) + 0.55 log(BAC)
GT = grazing rate (cells ml-1 hr-1) HNF = heterotrophic nanoflagellates (cells ml-1) T = temperature (C) BAC = prokaryote abundance (cells ml-1)
Start with non-zero grazing rate values, etc. Where… GT == Grazing effect hr HNF == Euk_ml BAC == PROK_ml T == highest temp
We only have these for MCR data.. can include GR for the “predicting euk_ml”
10.1 Predicted parameters
# head(all_vents)
allvent_wpredicted <- all_vents %>%
select(SITE, NAME, SAMPLE, EXP, PROK_ml, EUK_ml, GRAZING_EFFECT_hr, Highest.Temp) %>%
group_by(SITE, NAME, SAMPLE, EXP, PROK_ml, EUK_ml, GRAZING_EFFECT_hr) %>%
summarise(TEMP = mean(Highest.Temp)) %>%
distinct() %>%
mutate(PREDICTED_TEMP = (((log(GRAZING_EFFECT_hr)) + 3.21 - (0.99*log(EUK_ml)) - (0.55*log(PROK_ml))) / 0.028),
PREDICTED_GRAZING_EFFECT_hr = (10^(-3.21 + (0.99*log(EUK_ml)) + (0.028*TEMP) + (0.55*(log(PROK_ml))))),
PREDICTED_EUK = (10^(((log(GRAZING_EFFECT_hr)) + 3.21 - (0.55*log(PROK_ml)) - (0.028*TEMP))/0.99)),
PREDICTED_PROK = (10^(((log(GRAZING_EFFECT_hr)) + 3.21 - (0.99*log(EUK_ml)) - (0.028*TEMP))/0.55)))## `summarise()` has grouped output by 'SITE', 'NAME', 'SAMPLE', 'EXP', 'PROK_ml', 'EUK_ml'. You can override using the `.groups` argument.# View(allvent_wpredicted)sessionInfo()## R version 4.1.0 (2021-05-18)
## Platform: x86_64-apple-darwin13.4.0 (64-bit)
## Running under: macOS Mojave 10.14.6
##
## Matrix products: default
## BLAS/LAPACK: /Users/sarahhu/anaconda3/envs/r_4.1/lib/libopenblasp-r0.3.15.dylib
##
## locale:
## [1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8
##
## attached base packages:
## [1] stats graphics grDevices utils datasets methods base
##
## other attached packages:
## [1] gt_0.3.1 broom_0.7.9 cowplot_1.1.1 forcats_0.5.1
## [5] stringr_1.4.0 dplyr_1.0.7 purrr_0.3.4 readr_2.0.0
## [9] tidyr_1.1.3 tibble_3.1.3 ggplot2_3.3.5 tidyverse_1.3.1
##
## loaded via a namespace (and not attached):
## [1] Rcpp_1.0.7 lubridate_1.7.10 assertthat_0.2.1 digest_0.6.27
## [5] utf8_1.2.2 R6_2.5.0 cellranger_1.1.0 backports_1.2.1
## [9] reprex_2.0.1 evaluate_0.14 httr_1.4.2 highr_0.9
## [13] pillar_1.6.2 rlang_0.4.11 readxl_1.3.1 rstudioapi_0.13
## [17] jquerylib_0.1.4 checkmate_2.0.0 rmarkdown_2.9 labeling_0.4.2
## [21] munsell_0.5.0 compiler_4.1.0 modelr_0.1.8 xfun_0.24
## [25] pkgconfig_2.0.3 htmltools_0.5.2 tidyselect_1.1.1 fansi_0.5.0
## [29] crayon_1.4.1 tzdb_0.1.2 dbplyr_2.1.1 withr_2.4.2
## [33] grid_4.1.0 jsonlite_1.7.2 gtable_0.3.0 lifecycle_1.0.0
## [37] DBI_1.1.1 magrittr_2.0.1 scales_1.1.1 cli_3.0.1
## [41] stringi_1.7.4 farver_2.1.0 fs_1.5.0 xml2_1.3.2
## [45] bslib_0.3.0 ellipsis_0.3.2 generics_0.1.0 vctrs_0.3.8
## [49] RColorBrewer_1.1-2 tools_4.1.0 glue_1.4.2 hms_1.1.0
## [53] fastmap_1.1.0 yaml_2.2.1 colorspace_2.0-2 rvest_1.0.1
## [57] knitr_1.33 haven_2.4.3 sass_0.4.0